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SYNOPSIS 


Heat conduction problems with phase change have been 
studied extensively for over a century. These problems are 
encountered in a wide range of situations, encompassing 
such diverse applications as casting of metals and alloys, 
storage of thermal energy, spacecraft thermal design, drying 
of foodstuffs etc. In modern foundary technology, the 
prediction of the solidification rate and the temperature 
distribution during the solidification process is very 
important, in order to control the fundamental parameters, 
such as stripping time for static casting and with dr awl 
rates for continuous casting, as well as the incidence of 
casting defects. A series of analytical solutions have been 
reported in the literature, starting with the well known 
Neumann solution. Most of the investigators have so far, 
considered only conduction in the melt , neglecting the 
natural convection effects, which arise due to buoyancy. 

The study of natural convection in the melt is very important 
as it may affect the process in the following ways: 

(i) The motion in the melt may modify the heat 
transfer rate at the melt-solid interface, and hence affect 
the gross solidification rate, (ii) The convective motion 
in the melt may influence the structure of the solid formed, 
(iii) The distribution of the solutes in the multi- component 
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system may be affected by the bulk motion in the melt. 

(iv) The velocities in the melt -would affect the rate at 
which non-metallic inclusions can float to the surface. 

The present -work has been divided into three parts; 

(i) Computation of the solidification rate and 
the temperature and velocity distributions 
in a one dimensional solidification process. 

(ii) Computation of the flow pattern, solidification 
rate and the temperature and velocity distri- 
butions during solidification, in a two dimen- 
sional rectangular enclosure. 

(iii) Experimental verification of the analytical 
results obtained in (i) and (ii). 

In the one dimensional problem, equations have been 
written to study the solidification taking into account the 
heat transfer in the melt, in the solid formed and in the 
mold. Analysis has also been carried out to determine the 
velocity distribution that arises due to thermal gradients.' 

The equations are solved for different boundary conditions, 
such as constant wall temperature and convective and 
radiative heat transfer at the outer surface of tho mold. 
Solidification between infinite parallel plates end inside 
an infinite vertical cylinder have been considered. The 
temperature and the velocity distributions in the melt are 
obtained, along with the rate of movement of the solidification 
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front for various non-dimensional parameters such as the 
Prandtl number, Biot number, Stefan number etc. The computed 
results are compared 'with the results of other investigators 
for which comparison could be made. 

In the second part, the natural convection arising 
due to the thermal gradients during solidification in a 
rectangular enclosure hats been studied. The two dimensional 
equations have been solved, initially without considering 
the solidification as a chock on the numerical scheme 
and later, including the solidification. Bor the melt 
region equations for the conservation of mass, momentum and 
energy are written. These equations, coupled with the 
energy transfer equations for the solid region and the mold 
through the boundary conditions at the melt-solid and solid- 
mold interfere, are solved sinultanoously for different 
boundary conditions. The boundary conditions considered are: 

(i) The top and the bottom surfaces of the 

enclosure are adiabatic, whereas the outer 
surfaces of the mold are kept at a constant 
temperature . 

(ii) The top and the bottom surfaces of the 
enclosure are adiabatic while the outer 
surf ares of the mold exchange energy by 
convection and radiation. 

(iii) The top surface and the mold outer surfaces 
exchange energy by convection and radiation, 



the bottom surfo.ee of the enclosure 
being kept adiabatic. 

The non-linear, coupled equations are non-dimonsionaliz e& 
and solved numerically using the Alternating Direction Implicit 
Technique. The results arc computed for various dimensionless 
parameters, such as Rayleigh number, Prandtl number, Stefan 
number, aspect ratios, Biot number etc. Bor these boundary 
conditions, the temperature profiles, the velocity profiles 
and the rate of interface movement are obtained. 

Experiments were carried out, using paraffin wax as 
the solidifying material. The experimental arrangement 
consists of a rectangular Plexiglas tank with two aluminium 
plates inserted vertically, thus dividing it in three compart- 
ments. The molten paraffin wax is poured into the central 
compartment while the end compartments are used for circulating 
water at a constant temperature. The temperature of the 
paraffin wax is monitored by thermocouples, the output of 
which has been recorded, by a Hewlett Packard Data Acquisition 
System. The thermocouples are prepositioned, at given hori- 
zontal and vertical locations, from the temperature vs time 
plot, the position of the. interface at various time intervals 
is found. The experiments are conducted for different aspect 
ratios, and superheat of wax. The results are then compared 
with the theory. 
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i?T OI-iBtTG LATUR E 

Biot number at the mold outer surface = h-jd/k^; 

Biot number at the top surface of the enclosure - h^B/k^j 

specific heat [j/Kg °K] ; 

thickness of the mold [m]; 

diameter of non-metallic inclusions [m] ; 

•width of the enclosure [m] ; 

vertical velocity of the melt when Gr=l, [m/s]; 

heat transfer coefficient at the mold outer surface 
[j/s m 2 °K] ; 

heat transfer coefficient at the top surface of the 
enclosure [j/s m 2 °K] ; 

/*) 

gravitation constant [9«8 m/s^]; 

Grashof number = (gp (1^-2 C )Y^/ v 2 ); 

thermal conductivity [j/ s m °K] 

height of the enclosure [m] ; 

latent heat of solidification [j /Eg]; 

pressure [lT/m 2 ]; ZiydYOStcttic pressure C Wirt 3 - 

<d ynamte pressure c/v/zr^jy 
Prandtl number = v /a, 

heat transferred at the interface [j]; 

superheat present in the melt = p 0^ L Y o (T^-T c ) , [j] ; 

radial coordinate measured from 'die center of the 
cylinder [m] ; 

inner radius of the vertical cylinder [m] ; 

Rayleigh number = Gr.Pr; 

radiation constant at the mold outer surface 
= Ul cr d(T i -T c ) 3 )/k m ; 
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radiation constant at the top surface of the 

'Z 

enclosure = (£2° LCf^-ip^/k-jj 
Stefan number = L/C (1 .-1 ); 

p i c 

time [s]j 

Q O 

dimensionless time = vt/Y^ or at/D ; 
temperature [z] ; 

average temperature in the melt [Z] ; 
mold outer surface temperature or the ambient 
temperature as the case may be [Z] ; 
dimensionless temperature = 1 / (1 ^-T c ) > 
initial temperature [z] j 
temperature in the melt [iC] ; 
temperature in the mold [z] ; 
temperature in the solid [Z] ; 

dimensionless temperature in the melt = ( 1 c ) / i“^' c ) > 

dimensionless temperature in the moj-d =(V T c )/(T i" T c ) ’ 
dimensionless temperature in the solid=(T g -Tp/(f j 

solidification temperature [E] ; 
velocity in the x-direction [n/s] ; 

dimensionless velocity in the x direction = uD /oCjLf 
terminal rising velocitjr of noil-metallic inclusions 

Cm/ s] j 

velocity in the y direction [m/.s]; 

dimensionless velocity 1- in the y direction = vD/oc-^j 
vertical coordinate measured from the bottom of the 
enclosure Cm] ; 


x 
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x* dimensionless vortical coordinate = x/L; 

y horizontal coordinate measured from the center 

of the enclosure [m] ; 

y^ dimensionless horizontal coordinate in the melt = y/Y ! 

y m dimensionless horizontal coordinate in the mold 

= (y-T 0 )/(D-T 0 )} 

y dimensionless horizontal coordinate in the solid 

O 

= (y-i* )/ ); 

Y ! half width of the melt region [m] ; 

Y distance from the center lino to solid-mold inter- 

face [m] ; 

z coordinate measured from the center of the 

enclosure [m] ; 


Greek Symbols 

O 

a thermal diffusivity [m /s]; 

P temperature coefficient of cuibical expansion [K ; 

6 relative width of melt region = Y ! /Y ; 

relative width of mold region = D/Y : 
m o 

emissivity at the mold outer surface; 
emissivity at the top surface of the enclosure; 
p coefficient of viscosity [Kg/m s]; 

v kinematic viscosity = ji/p , [mg/s]; 


P 


p 

p 


P 

sat 


densit y [Kg/ m -3 ] ; 

•X 

densit y of non-mot allic inclusions [Eg/m -3 ]; 

density of the melt at the average temperature [Eg/m 
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a 

0 




* 


St of an~Boltzman constant [5,67x10""® J/s m^ K^] ; 
dimensionless t crap erature = ( 11 ~T ga ^_ )/ (T sa ^“T c ) > 
stream function [r / s] ; 
dimensionless stream function =4'D/a-^l; 
normalized stream function = $ / ^ ma2 -> 


4 > 


max 


Q 


opt 


to 


* 


maximum stream function [m®/s]; 
optimum relaxation factor; 
vorticity [s” 1 ] 

dimensionless vorticity = wD^/a^l; 
Superscript 

k iteration number; 

* dimensionless quantities; 

Subscript 

i grid point in x direction; 

j grid point in y direction; 

1 melt ; 

m mold ; 

n time; 

s solid; 

w boundary; 



CHAPTER 1 


IHTRODUCTIOIT 


Many processes involving phase change are encountered 
in chemical, food and metallurgical industries as well as 
in the fields of cryhiology and space technology# The 
applications of solidification or melting include, casting 
of metals and alloys, the production of chemicals as prills 
or flakes, the freezing and destruction of tumorous brain 
tissue, spacecraft thermal design and storage of thermal 
energy, etc. The prediction of solidification rate and 
temperature distribution during the solidification process 
is very important in modern foundary technology in order to 
control the fundamental parameters such as stripping timec 
for static casting and wit Mr awl rates of continuous casting, 
as well as the incidence of casting defects. Furthermore, 
solidification rate and temperature gradients, in the melt, 
the solid formed and the mold encompassing them, have a strap' 
influence on the final properties of the cast piece and or. 
the possibility of damage to the mold surface. Solidifies. cio 
rate and temperature gradients, however, are controlled by 
the direction and the rate at which heat is removed from 
the melt . 

The presence- of temperature variations in the melt 
will, in most cases, produce natural convection currents. 
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Under these conditions, the conduction problems in the 
solid and the mold have to he solved simultaneously with 
the natural convection flow and heat transfer in the melt. 
Historically, the possible existence of natural convection in 
the melt has largely gone untreated in the heat transfer 
literature. Therefore, many of the conduction-based phase 
change solutions developed by numerous investigators 
over a century do not fully reflect the physical reality. 

The inclusion of natural convection in the melt is very 
important as it affects the process of solidification in 
the foil owin g way s . 

(i) the bulk motion in the melt may modify 
the rate of heat transfer at the melt- 
solid interface aid hence affect the 
actual, gross solidification rate. 

(ii) convective motion may influence the 

structure’ of the solid formed, in particular 
the transit ion from columnar to equiaxed 
orientation, in case of solidification of 
metals, 

( iii) convection in the melt would affect the 

distribution of solutes in the multicomponent 
systems, 

(iv) velocities in the melt wotild affect the 
rate ot which non-metallic inclusions 
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can float to the surface, for ingot 
casting of metals. 

The objective of the present work is to study 
the thermal and fluid flow fields developed during the 
solidification, along with the movement of the melt-solid 
interface in an enclosed region. 

A brief review of literature on solidification 
and fluid flow in ai closed regions is given in Chapter 2 .' 

In Chapter 3 the temperature and velocity field developed 
during the solidification in one dimensional vertical slot 
and vertical cylinder are discussed. The dependence of 
solidification rate on various dimensionless parameters 
is also studied for the one dimensional solidification. In 
Chapter 4, the two dimensional solidification in rectangular 
enclosures with the natural convection in the melt has 
been discussed. Streamline and isotherm patterns are 
presented along with the position and shape of the interface 
at various times, for different boundary conditions. The 
experimental technique, apparatus and procedure are described 
in Chapter 5. The comparison of the experimental results 
with the theory is also presented in Chapter 5 . The 
conclusions and recommendations for the future’ work are 
given ih Chapter 6* 



CHA PTER 2 


LITERAT URE SURVEY 

An examination of the literature on phase 
change reveals no comprehensive or -unified treatment of 
the combined thermal and fluid flow phenomena that occurs 
during solidification. The phase change problems dis- 
cussed in the literatrire can be broadly divided into two 
categories; (i) Problems in which the conduction in the 
solid phase alone is considered assuming the melt to be 
at saturation temperature, (ii) Problems in which the 
superheat in the melt is considered. This superheat gives 
rise to "buoyancy forces, which in turn produce convective 
flows. 

2.1 COEDUCTION-BASED SOLUTIONS; 

2.1.1 Analy tical So lut ions ; Analytical and closed form 
solutions of the pure conduction problems with solidification 
have been obtained for a few idealized situations. The 
simplest case is the solidification of a serai- inf in it e medium 
with a constant surface temperature, known as the Stefan's 
problem. If the heat conducted into the solid alone is 
considered, the solidification time can be obtained by 
applying Chvorinov' s rule diseuissed by Cole (1969) • Carslaw 
and Jaeger (1959) have discussed the problem, when both 
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the solid and the melt are assumed to he seni— infinite, 
the melt being at a temperature different from the solidi- 
fication temperature. The solution to the above problem 
has been developed by Frans Neumann (1912) for the surface 
of the semi-infinite solid being maintained at a temperature 
lower than the solidification temperature. Tien (1965) 
applies the condition, that the surface temperature is 
time dependent, and suggests that the solidification thick- 
ness increases almost linearly with the time whereas the 
rate of solidification first increases and then decreases 
in a complicated manner. 

2,1.2 Honanalytical Solutions ; There are other methods which 
are less exact than the Neumann solution, and these 
generally require the use of a computer. A contour integral 
approach used by Boley (1964), graphical method suggested 
by Sarjant and Slack (1954), analog methods, applying 
electrical analogs to the heat flow, used by Faschkis (1945; , 
and isotherm migration method discussed by Dix and Cizek 
(1970) are some of the non-analyti cal methods reported in 
the literature. 

Another approximate but elegant method to solve- the 
conduction problems is proposed by Goodman (1958), employiuj, 
a heat balance integral technique. In the integral method 
treated by Hrycak (1963), the Neumann conditions are appl-Loo., 
but with an assumed temperature distribution which is linear 



6 


in the solid and parabolic in the melt. The heat flow 
at the surface of the solid, is removed from the system by 
convective flow to the ambient fluid. When the integral 
equation is specified, the equation can be solved using a 
computer. Hills and Moore (1967) have used the intogral 
method by assuming polynomial series for the temperature and 
the thickness of the solid formed. A generalized solution 
for the melting of a semi-inf inito solid, initially at the 
melting temperature and subjected to a non-linear boundary 
condition at the surface of the solid, has been discussed 
by Ozisik and Mody (1975), by applying the integral method. 

mile Goodman (1958) restricts the approximation 
of the temperature distribution inside the solid phase to 
a cubic polynomial, Ilegorlin (1968) extends the method to 
any ordered polynomial and suggests that a.s the order of rhe 
polynomial approximation is increased, the solution asympto- 
tically approaches the exact solution. Sparrow and 
Sham sunder (1974) have applied Mcgerlin’s technique to tho 
problem of the thermal energy storage in which tho heat 
conduction equation written for the cylindrical coordinate 
has been solved for tho convcctivo boundary condition. 

Tho other scheme which has created considerable 
interest in solving the heat conduction equation is tho 
numerical techniques. Goodling and Mhader (1974) have acpmic 
an implicit numerical scheme to solve the inward solidifi- 
cation problem for slabs, cylinders and spheres with a mixed 
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radiative and convective boundary condition at the outer 
surface of the solid. 


2.1.5 Treatmen t of the Hold s The previous solutions 
generally ignore the heat transfer in the mold and in most 
of the cases, the surface of the solid has always been 
assumed at some constant temperature. In real systems 
there is a mold that plays an important role in dissipating 
the heat to the environment . 

In a simple case, the mold is assumed to bo thick 
enough (semi- infinite) that no temperature rise occurs on 
its external surface, and initially, the mold is kept ax 
a particular temperature. In addition the somi- infinite 
molt is assumed to bo initially at its freezing temperature, 
thereby neglecting the superheat. A solution for the above 
mentioned problem has boon given by Schwartz (1935) and 
Ruddle (1957). 


2.1.4; Heat Transfer from the hold to Ambient : Heat transfer 
from the external surface of the mold determines the rate 
at which finite castings solidify and this has boon treated 
in several papers. Kosky (1975) has treated this problem 
using the concept of Biot modulus. Bn is the ratio of ol_- 
heat transfer coefficient at the surface of the mold to tx- 
thermal conductivity at the mold and is, thus, a measure ox 
the relative magnitude of the surface convection and internal 
conduction resistance to heat transfer. (A low Biot number 
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implies no internal conduction resistance, so that the 
temperature is uniform throughout the solid,) Szekoly and 
Ludely (1966) have considered the surface flux as a mixture 
of convection and radiation from the surface of the solidi- 
fying molt. An art itrary flux at the surface is treated 
by Hamil and Bankoff (1963), using the integral method. 

G-arcia and Prates (1978) have developed a model for the 
unidirectional solidification of the metals in the molds 
cooled by fluids such as air and water, Their model permits 
the measurement of the Newtonian heat transfer coefficient 
at the metal/ mold interface, end it differs from other 
analytical models in the sense that it is more general in 
application and easier to mcnipulc..to, while retaining the 
advantage and convenience over numorica.1 techniques, 

2.2 THERMAL CONVECTION; 

The movement of a fluid in free convection results 
from buoyancy forces imposed oil the fluid, when its density is 
increased locally as a result of the cooling process. The 
two main sources of the convective flow, that can be 
considered are: (i) the temperature difference between the 
hot fluid and some colder wall, either the melt— solid 
interface or, more generally any cold vortical wall immersed 
in the hot fluid, (ii) the heat loss at the top free surf .ace 
of the fluid exposed to the environment. Radiation and 
convection to the ambient will locally cool the melt near 



tho top suxfn.ee, and. the heavier fluid then moves downward, 
displacing the warmer fluid "below. 


2.2.1 Convection without Solidification : Numerous investi- 
gators have considered tho presence of natural convection 
of the fluid, in an enclosed region hut without solidification. 
A study of those investigations is helpful to understand tho 
flow pattern that nay arise with solidification. Eckert 
and Carlson (1961) have provided very useful information on 
tho temperature profiles and the circulation pattern developed 
hy an air layer, enclosed between two isothermal vertical 
plates. Elder (1965) has analysed the temperature and 
tho velocity profiles of tho fluid enclosed in a. narrow slot 


and has arrived at tho criteria. 


tho establishment 


secondary and tertiary flows. The application of tho finite- 
difference numerical techniques has boon pioneered by 
Wilkes and Churchill (1966) who have studied tho steady state 
and transient behaviour of the fluids enclosed in rectangular 
cavities for various boundary conditions. In a subsequent 
work Torrence and Rockett (1968) have made tho transient 
analysis for the fluid in a rectangular enclosure heated from 
below, for high Rayleigh numbers, whereas tho work by 
Elder (1966), Do Vahl Davis (1968) and Nowell and Schmidt 
(1970) has concerned with the steady state systems. 


2.2,2 Solidification with Conv ection; Although the natural 
convection heat transfer and solidification have been 
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extensively studied individually, their combined occurrence 
has so far received little attention* Tiller (1962) has 
recognised the role of convection in the early stage of 
solidification. Ostrach and Kroger (1974) have considered the 
role of natural convection in the analysis of the steady 
state phase change problem associated with continuous casting. 
Seki, Pukusako and Sugawara (1977) have investigated the 
effects of the natural convection in the problem of freezing 
or melting of water in horizontal layers, using linear 
p erturbat i on t e ehniqu e . 

While the works mentioned above represent the effects 
of fluid motion, the first attempt to study the exact nature 
of the fluid motion has been made by Szekely and Stanek (1970). 
They have developed explicit expressions for the temperature 
and the velocity distribution in the fluid. However, since 
the relevant energy equation does not contain any velocity 
terms, the heat transfer process is not affected by buoyancy 
and conduction is the only transport mechanism. The only 
investigation which takes into account the role of the 
natural convection in the energy transport is that of Sparrow, 
Patankar and Ramadhyani (1977). They have carried out an 
analysis for the two dimensional melting in a cylindrical 
enclosure, taking into account the natural convection induced 
by temperature differences in the melt. Their results 
decisively differ from 'those of the earlier conduction 
solutions. 



•LJL 

2.3 ESPERIHH\PTSs 

In addition to the above mentioned investigations, 
there appears to he a greater awareness in other literature, 
particularly metallurgy, ox the role of natural convection 
in phase change problems. Experiments have been conducted 
by many investigators to understand the role of the natural 
convection in modifying the structure of the solid formed, 
and distribution of the solute in the multicomponent system. 

In an experimental study Cole and Bolling (1965) have shown 
that the structure of the casting can be significantly 
changed by eliminating the natural convection effects. An 
attempt to know the quantitative effect of the natural 
convection has been made by Szekely and Chhabra (1970) who 
studied the controlled solidification of lead under conditions 
of unidirectional heat flow within the system. The observed 
shape and position of the melt-solid interface agreed reasonably 
with the predicted results after giving consideration to 
the heat transferred by natural convection from the melt to 
the solid. The experiments have been performed for both 
steady state and transient cases. Sparrow and Ramsey (1978) 
have conducted experiments to provide definitive information 
about the heat transfer processes which occur when the solid 
is melted by a vertical, cylindrical heat source embedded 
within the solid. The role of natural convection has been 
studied by talcing photographs of the melt-solid interface 
at various time intervals. Recently Sparrow, Ramsey and 
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* Kemink (1979)" have performed experiments for freezing 
n-eicosane, a paraffin, under conditions where the melt is 
either above or at the fusion temperature, the freezing 
material being housed in a cylindrical containment vessel. 
Interesting results have been obtained when the melt is 
superheated, for the boundary conditions considered by 
the authors. The other experiments intended to clarify the 
interaction of the natural convection and the phase change 
in one dimensional geometries are reported by Heertjes, 
Jongenelen, and de Leeuw den Bouicr (1970), and Chiesa 
and Guthrie (1974). 



CHAPTER 5 


ANALYSIS OP SOLIDIFICATION IN ORE DMENSIONAL 
VERTICAL SLO T AND CIRCULAR CYLINDER 

3.1 INTRODUCTION: 

The mathematical modelling of the solidification 
of the nelt, in molds of finite dimensions is very complesc. 

The solution of the equations, to understand the thermal 
and fluid flow phenomena, is complicated since in a multi- 
dimensional system the thermal and momentum equations are 
nonlinear and coupled, and have to he solved simultaneously 
along with the moving boundary equation, as the solidifi- 
cation proceeds. Solidification in one dimensional, vertical, 
narrow slot (or vertical cylinder) represents the realistic 
situation, when the height of the mold is very large 
compared to the width, i.e., molds of very high aspect 
ratios. Because of the one dimensional approximation, the 
thermal and the momentum equations are decoupled and 
conduction is the only mode of heat transfer. Even for such 
a simplified system analytical closed form solution is not 
possible as the hc-at conduction equations for the melt, 
solid and the mold have to be solved simultaneously along- 
with the moving boundary equation. 

In the present chap ter, an analysis is carried out to 
study the transient thermal characteristics of solidification, 
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accompanied by natural convection flow, in a vertical 
slot and in a vortical cylinder. The solution is facilitated 
by suitable coordinate transf ormat ions that immobilise the 
moving interface between the melt and tho solid. Consideration 
has been given to the energy conducted through and stored 
in the mold, solid and the melt regions. The problem has 
been studied for various boundary conditions at the mold 
surface. These are; 

[i] constant wall temperature; 

[ii] convection from tho mold to ambient, and 

[iii] radiation from the mold to ambient. 

The governing coupled equations for the melt, solid 
and the mold regions are solved rising an implicit, finite- 
difference scheme to obtain the time dependent temperature 
and velocity distributions. The effect of various non- 
dimensional parameters such as Prandtl number, Stefan number, 
Biot number etc., on the rate of solidification has also 
boon studied. It has been shown that in case of liquid 
metals, tho velocity of the melt that arises due to thermal 
gradients is sufficient to float up any non-metallic 
inclusion of dimensions found in practical applications. 

Though the analysis has been made for the parallel 
plates and the circular cylinder, only tho parallel plates 
equations are discussed in detail; however a brief idea 
about tho governing equal ions for the circular cylinder and 
their differences from tho equations of the parallel plates 
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have boon discuss od in section 3.4. 

3.2 ANALYSIS s 

A schematic diagram of the vortical parallel plates 
and the circular cylinder is shown in Figure 3.1. The 
width of the rectangular slot is 2 Y q and the inner diameter 
of the cylinder is 2B. Q . The slot is extended to infinity 
in both the x and z directions, whereas axisymmetry is 
assumed for the circular cylinder which is extended to 
infinity in the x-diroction. The thickness of 'the mold is d. 

At time t=0 , the mold is kept at a temperature T , the 
ambient temperature. At this time, melt at on initial 
temperature T^ is poured into the mold. As a result, thermal 
gradients arc set up and a thin crust of solid is formed at 
the inner surface of the mold. An inward movement of the 
melt-solid interface also starts. The transient temperature 
profile developed in the melt changes the local density of 
the melt giving rise to buoyancy forces, which in turn cause 
natural convection in the melt. 

3.2.1 Assumptions ; The following assumptions are made for 
the analysis: 

[i] The physical properties such as thermo,! conduct ivitj", 
specific heat and density are o.ssujned to be independent of 
the temperature for the molt, the solid and the mold with 
the exception of the density of the melt contributing to the 
buoyancy forces. This assumption is well justified, as 
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long no the temperature? variation is not high. 

[ii] The melt is considered to he Rewtonian. This 
assumption is also valid as in most of the practical 
applications, the fluid is Rewtonian. 

[iii] The vortical slot is considered to he narrov 
so that the heat transfer can occur only in the y-diroction. 

An immediate consequence of this assumption is that only 
the vertical velocity component u predominates, leaving the 
horizontal velocity component v negligible, since parallel 
flow arises for this configuration. The problem becomes 
unidirectional and hence the solid front con move only in the 
y— direction, parallel to the mold. 

Bo cause of this assumption the- momentum and the 
energy equations one decoupled md conduction becomes the 
only mode of heat transfer in cooling the molt. Batchelor 
(1954) has discussed that the heat transferred between two 
fixed, vertical plates can bo assumed to bo by conduction 
provid od 

Ra <500 (5.2.1) 

where Ra is the Rayleigh number. Though the above relation 
is very approximate, it gives an idea about the mode of heat 
transfer. 

[iv] The convective flow developed in the melt is 
assumed to be laminar. The criterion for the onset of 
turbulent natural convection has boon discussed by Colo (1969) 
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who puts the laminar range as 0 < Ra < 10°; the assumption 
is quite vaUd since the present investigation is well 
within the range. 

[v] Time taken to pour the melt into the mold is 
small and neglected in calculations, so that solidification 
starts only after the pouring is complete. 

[vi] Solidification is assumed to occur at a 
constant temperature. 


3.2.2 Governing Eque,tions; The problem can be stated 
mathematically, by the use of the following equations. 

(i) Equation of continuity in the molt. 

(ii) Momentum equation in the molt . 

(iii) Energy equations in the molt, the solid and 
the mold. 

The equations aro written for only one half of the 
slot a.s they aro valid for the other half also, due to 
symmetry. 

3.2.3 Energy Equations; The energy equations arcs 


Melt ; 


3 T. 
Tt* 


,2 T 


= a- 


sy 


St 


= a 


a 2 g s 
9y 2 


(3.2.2a) 


Solid : 


(3.2.2b) 
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Mb Id: 


3 1' 


JSL 


>2 1 


a 


9 1 “ “m " ; 2- 

9y 




(3.2.2c) 


whore 1^, l g and T are the temperatures of the no It, 13a e 
solid and the mold respectively. 

The initial and the boundary conditions for these 
equations can he written as: 
in th e nolt , 


t = 0: 


y = 0: 

y = T« : 
where 1 


1 


i 


i 


(3.2,3a) 


-= 0 (because of symmetry) (3.2.3b) 

(3.2.3c) 


nd, 


sat 
3 1 


3 y 

T — 1 
X 1 “ x sat 

= solidification temneraturo 


31 . 


k 


s 3 y 


1 37 


pi 


dY* 
QL dt , 


(3.2.3d) 


1 = l 

f3 Sal u 


(3.2.3o) 


in the solid , 

y = I* s 

y = x : 

J o 

At the nolt— solid interface the heat relee,sed by the 
solid should bo conducted through the mold. Hence, 
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s 

in the mold . 

t = 0: 

y = Y : 
J o 

y = I: 


~“® ' “ k 
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i 


n 


T 

< 

1 


(3.2.3f) 

(3 . 2.3g) 
(3.2.3h) 
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Boundary condition of first kind, 
= T q = c one trait 
(or) 

Boundary condition of third kind. 


3T 

-k — Ja = h (2 - gj ) 

ei 5 y 1 V El c' 

(or) 

roundary condition of fourth kind, 

3 T y 
_1 

n ay 


-k 


s = e^T 


T 4 ) 

G 


(3-2.31) 


3*2.4 I'l o.nc£. it nn B guati on: 

Since the horizontal velocity v(:s,t) is negligible 
in couparioon with the vertical velocity u(y,t), the 
momentum equation for the u-comp oneiit can bo written as % 


3 u 3^u 

p Tt = ^ . 

3y 


2 ~ PS. 


.±£ 

3 r 


(3.2.4) 


Tlio firct torn on the PUS ic duo to momentum diffusion 
whereas tho second raid the third tome arc duo to body forces 
and tho pressure respectively. 

Tho second raid the third tons can be comb in od into 
a single term, tho details of which are discussed in 
App ondix I . Thus , 


3 u 


3 t 


v 2 + Pg(T-T) 

ay 


(3.2.5) 


The temperature T ic an average temperature in tho molt, 

the value of -which is found in such a way that tlio cons orvat ion 
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of mass is satisfied. 

The initial and the boundary conditions for this 
equation can bo written as: 



■u = 0, 

at 

t = 0 

(3.2.6a) 


—is = 0 , 

9 y ’ 

at 

y = o 

(3.2.6b) 

and 

u = 0 

at 

II 

(3*2. 6c) 


3*2*5 Continuity Equation : The continuity equation for 
the mult is written as: 



(3.2.7) 


since the velocity components v and w arc neglected in 
comp arisen with the u~conponm.it velocity. Apart from the 
assumption (i) if wu also assume that there is no change* of 
volume during solidification then the above equation becomes, 


— - = 0 ( 3 . 2 . 8 ) 

9 x 

Thus, the density of the melt is assumed to bo 
constant in the continuity equation, whereas it is assumed 
to vary with the tcr.ipora.ture in the equation of notion. 

This kind of treatment is very common in natural convection 
studios and it has boon discussed by Schlichting (1969). 

¥ith those simplifying assumptions, which arc thought 
to be reasonable at least as a first approximation, the 
equation of continuity nay bo stated by expressing the fact 
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that tiior o is no not notion of fluid across any horizontal 
section, i t o. 

Y* 

/ udy = 0 (3.2.9) 

o 

Equation (3.2.9) is a nodifiod fom of continuity 
equation, which can "bo used directly in the computation 
and has t><--en discussed in the section 3.3.2. 

3.2.6 I Jon-dimensi onalization ; Since the analytical closed 
form solution is not possible, the problem is solved by 
adopting numerical techniques. The problem can be solved 
with less complication, if the moving boundary is immobilized. 
To facilitate this, a transformation of the variables is 
carried out. The other purpose of the transformation is 
to introduce dimensionless variables and parameters. 

The non-dimensional variables for "the melt are: 

y* = y/Y* (t ) 

T* = (Ti-T c )/(T i -T c ) (3.2.10) 

u* = u Y q /v 
E or the solid, 

y* = (}--Y'{t))/(Y 0 -Y'(t)) 

h = - T c )/(T i - T o } 

and for the mold. 



(3.2.11) 
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4= (y-V/ (D - T o> 


T* = (T 

t-i ' i 


n -n ‘ V^ T 1 " V 


(3.2.12) 


lho dinonsionlcqs tine is defined as, 


t* = vt/Y 2 
7 0 


(3.2.13) 


Because of the t ran of o m •. ax ion, flic dinonsions for 
the domino of the nolt, the solid and tlio uold bocono 
•unity, and do not change with tine. This sinplification, 
however, results in greater coup lenity in the governing 
equations ,ns shown below. 

By substituting these dinuns ionic ss variables in 
the- governing equations (3.2.2a) - (3.2.2c) wo got, 


Melt : 
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3T 


1 


. _ r * 
3 /i 


y* 
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do ... 

dt* 


1 

Pr c 2 


6 (t 4 ) 


2 T* 


s y' 


*2 


(3.2.14a) 


where. 


6(-fc*J = Y‘(t)/Y 


o 


and. 


Pr = Pra.ndtl number 


Solid: 


3 t ' 


* 


aT * y -1 


3 7, 


* (i-o Ct^T 


d6 

dt* 


a 2 

° S 


1 pr (l-dtt 4 )) 2 3 y* 2 


(3.2.14b) 


and, 
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Mold: 



(3.2.14c) 


whore, 

6 = D/Y 

n 1 o 

Thu extra tor no found in tim LHS of equations (3.2.14c) 
mid (3. 2.14h) represent the convection associated with the 
innob iliz at ion of the noving boundary. An observer sitting 
on the ir mobilized nolt-solid interface sees a mss noving 
towards (or away) fron him 

The boundary end the i n itial conditions b econo: 


in the nolt 


*i- 1 

3 T* 

ayj 

m * _ _ 1 

X 1 “ 1 + 0 


at t*= 0 

at = 0 

at 3 ^ = 1 


whore, 0 = (It - T 


sat 


)/(T sat -T c ) 


(3.2.15a) 

(3.2.15b) 

(3.2.15c) 


= Constant 



6TtT 




Pr.Stc 



(3.2.15d) 




(3.2.15o) 

(3.2.15f) 

(3.2.15g) 

(3.2.151a) 




26 


and T* = T c /(T ± - T q ) 
(or) 


1* = 0 

n 


for the ■boundary condition of first kind 


(3-.2.151) 


Thu nonontun equation gets transformed to 
.* 


3 V 

vrh ere , 


3u* _ au* ^1 d6 _ _1 

„, r * FpyJ df* - 2 


o (t*) 3y 


«2 * * 

2—3 + &r (T* 


*2 


• n 
(3.2.16) 


Gr 


g P (T- 


1 ) J D 

c ' o 


= Grashof nunbe: 


The initial and the boundary conditions becone, 

u*= 0 at t*= 0 (3.2.17a) 

= 0 at y * = 0 (3.2.17b) 

3 7i 

u*= 0, at y * = 1 (3.2.17c) 

By applying the sane non-dimensional parameters, the 
continuity equation (3.2.9) can be transformed to 
1 

/ u* ay* = 0 (3.2.18) 

o 


3.3 SOlUTIOtT IIETEODOLOGTs 

The equations (3.2.14) - (3.2.18) have been solved by 
numerical techniques. As a first step, the thickness of the 
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solid formed, (1-6), was obtained from the moving boundary 
equation (3.2.15d). Equation (3.2.15d), when written in 
finite-difference form, reduced to cubic polynomial, which, 
was solved using Ch.eb 2 rsh.ev approximation technique discussed 
by Carnahan, Luthur and ¥ilkes (IS 69). 

At time if = 0, the starting difficulties were 
avoided by assuming a very snail thickness of the solid, 
in this case taken as 0.0001 Y . The temperature profile in 
the solid was assumed to be linear at time t* = 0, which has 
been supported b 2 r Schlichting (1969) . Verification of this 
assumption was also made with Neumann’s solution, which 
one can expect to be valid at the initial stages. The 
initial solid thickness calculated by the present scheme 
was compared with the Neumann’s solution and agreement ~o 
within 0.1 per cent xtfas obtained. The initial thickness of 
the solid laxrer, due to chilling, was also varied to ascertain 
that it does not affect the numerical results significantly. 

The second step involved in the computation is the 
simultaneous solution of Eqs (3.2.14a) ~ (3. 2. 14c)* The 
equations wore solved using the finite-difference technique. 

The melt, the solid, and the mold regions wore divided 
separately into 10 grid spacings each, along the y-axis, 
each grid being of length 0.1. The partial derivatives of 
the equations were replaced by their finite-difference 
approximations , which convoru 'ohe governing partial differential 
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equation 3 into a sot of simultaneous algebraic equations, 

3.3*1 Explicit Scheme ; The finite-difference explicit 
scheme, discussed by Carnahan, Luther and Wilkes (1963) 
when tried, was found to be unstable, The instability 
is duo to very small domain in the solid region, in the 
beginning of the solidification, which has restricted the 
time step to an unreasonably small value, 

5*3.2 Imp li c it Scheme : The implicit scheme, which is 
found to be suitable for this problem, and more accurate 
than the explicit scheme, overcomes the difficulties 
encountered in the explicit scheme, at the cost of a some- 
what more complicated calcuiational procedure. It consists 
of representing a derivative 3 1^/ 3y£ by a finite - 
difference form, evaluated at the advance point of time 

t* , t instead of t* as in the explicit scheme,, 
li+l n 

As an example, Eqn. (3, 2,14a) can be written as 



(3.3.1) 


This can be simplified to, 





* * 


X , l T l,j-l,n + l + ( 1+2 * V l + E 1 1 ” > T l,j,n + 1 


-( x ;j_ Ay* + x q) t i, 3 + i iI1+ i - T i,j, n (3.^-2) 


where X — At*/( Ay*) 
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5 n+l 


Equation ( 35 , 3 * 2 ) represents a set of simultaneous , 
algebraic equations fori! in g a tridiagonal matrix, 

The tridiagonal matrices, thus formed for the melt, 
the solid and the mold, were first solved by the conventional 
Gauss elimination method. One of .the major disadvantages of 
the Gauss elimination method is the accumulation of the round- 
off errors, especially- at higher time steps. In view of this, 
this method was only used to solve the equation of the melt 
region. Eor the solid and the. mold regions, the governing 
equations along with their linear boundary conditions were 
written in 'a matrix form and the matrix was inverted using the 
maximum pivot strategy. In case of the radiative boundary 
condit i on ( Eon ,(3.2.15i)) whi eh is non— linear , simp lif i cat i on 
has been made to bring it to a linear form, as discussed 
in Appendix II. 

Interestingly, though the computational time involved 
in matrix inversion using maximum pivot strategy is much more 
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than solving the tridiagonal 
technique, the former allows 
wich reasonably good accuracy 


matrix using Gauss elimination 
the use of higher time step 

t 

* Hence, the computation time 


lost 


in the matrix invers 


on can he compensated by using 


large time steps. The At* used in the computation is 0,01, 
in wiij-Cii case a CPU xime of 1 minute is needed for 70 per 


cent solidification on a DEC 1090 computer. 


The momentum equation (3,2,16) which is similar to 
the energy equation (3,2,14a) has "been solved Toy the Gauss 
elimination technique. The important point here is that 
the quantity T * in Eqn.(3.2.l6) is fictitious and has to he 
adjusted iteratively for each time stop, in such a way that 
the continuity equation (3.2,18) is satisfied. 


3.4 VERTICAL CYlIiTDEHs 

The analysis has also been done for the vertical 
cylinder. The governing equations ares 
Energy equation; 




9 

Tr 


(r 


9 T, 


? r 


where k = 1, for the melt 

s, for ‘the solid 

m, for the mold 

Momentum equation; 


(3.4.1) 


3 U 

at 



3 

9 r 


(r—g) + gi3(T-T) 


(3.4.2) 
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The as sump t ions, boundary conditions, non-dimensionali- 
zation and the method of solving these equations are the 
same as those for the parallel plate case, except for the 
slight change in the governing equation, due to the intro- 
duction of an additional tern in the r-coordinate. 

3.5 RESULTS MTD DISCUSSIOITs 

The presentation of the results is mainly from two 
view points. 

(i) The study of the temperature and velocity 
profiles. 

(ii) Dependence of the interface movement with time 
on various non-dimensional parameters. 

3.5.1 Te mperature Profile : figures 3.2 and 3.3 indicate the 
t omperatiire profiles obtained for the melt, the solid and 
the mold regions, with dimensionless time as parameter, for 
the constant wall temperature and the convective boundary 
condition respectively. The curves lying above the solidi- 
fication temperature represent the melt region whereas those 
lying below tho solidification temperature, hut before tho 
melt-solid interface represent the solid region. This type 
of representation is logical as it gives an idea ahotxt tho 
position of the melt-solid interface at a particular time. 

The dotted linos shown in figure 3.2, indicate the temperature 
profiles for the cylindrical case. In the figures 3.2 and 
3.3, for the non-dimensional parameters considered, the 
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temperature profiles for the solid and the mold are linear. 

This is due to very high thermal diffusivities of tlio solid 
and the mold, which lead to insignificant storage terms in 
their energy equations. It is seen from figure 3.2 that 
the temperature decays at a fastei: rate in the cylindrical 
case than that in the planer case, resulting in a more rapid 
solidification. This fact is hotter observed from Figure 3.5, 
where the thickness of 'die solid formed can be compared for 
the planer case with that for the cylindrical case. 

3*5*2 Velocity Profile ; The time dependent velocity distri- 
bution in the melt is plotted in Figure 3.4. A point which 
c an be immediately understood from Figure 3.4 is that, certa,in 
amount of time is required for the development of the 
velocity field. At a given value of y/D, the maximum velocity 
is reached at some intermediate -erne. At larger trmos as 
the solidification proceeds, the velocity field decreases. 

An examination of Figures 3.2 end 3.4 indicates that 
for low Prandtl number the fluid motion decays more slowly 
as compared to the thermo 1 field. This suggests, that 
apu re ci able fluid motion con persist even after mosc of one 
superheat has been removed from the system. This behaviour 
is inherently associated wrth the low Prandtl number fluids 
(e.g. liquid metals), with the physical fact that the moluen 
metals con transmit heat much faster than they are able to 


t r an sf or m om ontum 



Vert ica l 
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The velocity profiles are symmetric about y—0 and 
the total upward flow equals the total downward flow, so 
xhat an enclosed region is simulated# 

Inspection of Eqn.(3.2.16) suggests that the nelt 
velocity varies with the Grashof number as, 

u*(y,* t*) = Gr. f(yj t*) (3.5.1) 

where, f(y* t*) = velocity of the nelt when 
Gr = 1 

This has been verified by computing the velocity 
u*(y* t*) for different Grashof numbers. Thus one can 
expect higher velocities at higher Grashof numbers. The 
width of the slot X Q , plays an important role in the magnitude 
of Grashof number as it is raised to its third power in the 
definition. This suggests, in the case of larger slots, 
the natural convection developed can be significant even if 
the superheat of the nelt is low. 

3*5.3 Effect of Melt Velocity : The effect of nelt velocity 
can be understood by calculating the distance travelled 
inclusions in the nelt. Figure 3.4 represents the velocity 
distribution for the typical case of ingot casting of steel 
(properties in Appendix III) in a slot of 1.0 cn width. 
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Pro n the figure, the maximum velocity of the melt u is 

found to be 5.35, which is equal to 1.02 x 10 ^ n/sec in 

absolute values. In this case a superheat of 400°C is 

assumed which is hypothetical. In actual case the superheat 
be 

will/ about 20 C and the nold width is In, However, for 

the actual case the velocity u* is much nore than 5.35 as 

8 

the Rayleigh number is 6 x 10. The calculations are not shown 
for this case as it runs contrary to the assunption(iii) of 


Sec, 3. 2.1. , . 

If we assure that the molten steel contains non- 

metallic inclusions of average diameter 2.0 microns then 

tlio terminal rising velocity is calculated by the Stoke' s 

law; 


rz 

■ ( p - P ) g d/ 

- 18 n 


(3.5.2) 


where, 

3 3 

Pp = density of the particle =3.0 x 10 kg/n 

p = density of the melt 

d — diameter of the particle 
P 

Thus , 

C 

= 1.4 x 10“ n/sec 

This value is much smaller than the maximum velocitjr 
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of the no It , which, implies that the particles would move 
along with the buoyancy induced flow that arises duo to 
the temperature gradient in the melt. 


3-5.4 E ff ect of Prantl H umber (?r).; The thickness of 
the solid formed with respect to time, with Prandtl number 
as the parameter is shown in figure 3.5. The rate of solidi- 
fication is found to increase with a decrease in the Prandtl 
number. This is expected as fluids with a low Prandtl number, 
have a higher thermal diffusion, resulting in large heat 
transfer and hence a faster solidification. It is seen from 
figure 3.5 that the rate of solidification is linear at larger 
times. At large times, most of the superheat has been removed 
from the melt and the molt temperature reaches approximately 
the solidification temperature. In such circumstances only 
the latent heat is to be removed. Hence the rate of beat lib era o eel 
at the melt-solid interface is approximately constant. The 
heat extracted through the outer surface of the mold to the 
ambient also remains constant as the storage terms in heat cono.uc c 
equation 4 s vecy small. This can bo sc-en from figure 3.3. 
which shows the temperature gradient to be almost constant wiidi 
time. Since the heat liberated at the interface as well as the 
heat removed at the mold surface remain constant, the interface 
growth rate is linear at larger times. 

3.5.5 Effect of St of an Humber ( SteJ_ : The effect of Stefan 
number on solidification is shorn in figure 3.6, in terms of 
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•6 Dependence of interface movement on 
Stefan number for the boundary 
condition of third kind for vertical parallel 
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the dependence of the thickness of the solid formed on 
tine elapsed. It is found that the solidification proceeds 
at a faster rate when Sto is snail. This is true as small 
Sto implies little latent heat and hence a faster solidifi- 
cation rate. 

3.5.6 Effect of Superheat (0) ; The movement of the- solid 
front with time, with 0 as a parameter is shewn in figure 3.7* 
It is seen from Figure 3.7, that the rate of movement of the 
solid front decreases with increase in 0 values. This is 
reasonable as physically a high superheat (or largo 0) 
implies more boat load and hence the decrease in solidifi- 
cation rate. 


3*5*7 E ffect of Biot number, Rad iation Constants; The 
importance of Eiot number Ei^ end the radiation constant 
p t c 1 can bo understood from figure 3.7 and 3.8 respectively. 
The thickness of the s olid formed at various times is found 
to increase with Bi 1 (or llc^ , as a high value of Bi 1 (or Rc 1 ) 
would imply a high convective (or radiative cooling), leading 
*fco rapid solidification as oTosorvoa# 


3.5.8 C ompariso n with Previous Work : figure 3*9 shows the 
results computed by Szokcly and Stanok (1970), IToumann (1912) 
and the present work (Ramachandron, Jaluria and Gupta (1981a)), 
for a 'topical case of solidification of stool. The solution 
has boon computed neglecting the resistance offered by the 
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mold, The* results are in good agreement with Szekely and 
St an ok (1970) who assumed a linear temperature distribution 
in the solid and neglected the effect of mold. The solution 
predicted by Neumann shows a decrease in s olid if i cat ion 
rate as Neumann assumed semi-infinite behaviour in the solid 
and in the melt. 


3.6 COTTCLUDIUGr REMARKS : 

An analysis has boon carried out for the beat 
conduction problem during solidification between parallel 
plates and in a vertical cylinder. The heat conducted through 
the melt, the solid and tho mold as well as the velocity 
field generated in the molt duo to thermal gradients, have 
been considered. Tho effects of various non-dimensional 
parameters on the solidification rate have also boon studied. 

The one dimensional analysis is valid only when tho 
height of the mold is very large- compared to its width. 


Chapter 4 describes analysis for 
fi cat ion, for which the height to 


the two dimensional. solidi- 
tho width ratio is small. 



CHAPTER 4 


.yWO.PjnMTSIOirAL SOIIDIFICATION INSIDE a rtEC-l'AITGULAR ENCLOSURE 


4.1 INTx. ODUC 1 1 ON : 


Hie solidification ox melt in enclosures of finite 
dimensions lias not "beer extensively studied in the literature, 
since the mathematical formulations and their solutions are 
extremely difficult and 'many simplifying assumptions have 
to he made. Hany of the solutions available now are 
conduction-based and they do not physically reflect the 
reality. The most important terms in the equations which 
define the tt.ro dimensional solidification problem are the 
convective terms which naLe the equations non-linear and 
rule out the possibility of any analytical solutions. Even 
to solve the equations numerically, special formulations are 
needed as non-linear derivatives lead the problem to a. 

_ '■ . numerical instability. In most of the cases 

the convective terms cannot be neglected as internal thermal 
f loirs of varying degrees are almost always present in the 
melt, since even slight temperature differences in the melt 
can give rise to buoyancy forces which in turn can produce 
sizable convective flow. The flow thus developed is important 
in determining the shape of the interface, solidification 
time and the solute distribution. Further the accounting 
of natural convection exposes the differences between solidi- 
fication and melting problems whose conduction -based solutions 



are interchangeable. Thus the problem involves the 
solution of nonlinear, coupled, simultaneous equations 
wrioten for the thermal and flow field in the melt and the 
energy equations written for the solid and the mold. 

In the present chapter an analysis is made to 
under s tana quantitatively the thermal and flow field 
developed during solidification in a rectangular enclosure, 

Two dimensional equations describing the temperature and 
flow field are written. Consideration has also been -given 
to the heat conducted and stored in the solid and the mold. 

The rate of interface movement which is unknown a priori 
is calculated by solving the heat balance equation at the 
molt— solid interface. The problem is solved for various 
boundary conditions. 

The coupled partial differential equations arc non- 
dimcnsionaliscd and solved numerically using Alternating 
Direction Implicit technique and Successive Over “Relaxation 
method. The equations are solved, initially without consi- 
dering solidification, as a check on the numerical scheme 
and later including solidification. The results arc computed 
for various dimensionless parameters such as Rayleigh number, 
Prandtl number, Stefan number, aspect ratio etc. The signi- 
ficance as well as the effect of these parameters on solidi- 
fication arc discussed. The temperature and velocity profiles 
arc; also obtained along with the rate of interface movement. 



4.2 FORi-iUIATIOITs 

Consider a rcctanG^-lar enclosure of width 2Y Q , 
mold thickness d and height 1. The enclosure is extended 
to infinity in the z-diroetion as shown in Figure 4.1. 

Initially, the nold is kept at a temperature T which 

V*/ 

is helov the solidification temperature 1 , At tine t=0, 

S cl x • 

the molt at temperature T. ( > T .) is poured into the 
enclosure. Immediately, a thin crust of solid is formed at 
the inner surface of the mold and an inward movement of 
solid front starts. Because of the thermal gradients set 
up in the melt, convection currents start in the melt. 

4.2.1 B oundary C on ditions ; The following cases are consi- 
dered: 

CASE I: The top and Bottom surfaces of the enclosure 

are kept adiabatic whereas the outer surfaces 
of the mold are kept at a constant temperature. 

CASE II: The top and bottom surfaces of the enclosure 

are kept adiabatic while the outer surfaces of 
the mold exchange energy by convection and 
radiation with the surroundings which is at a 
t emp eratu.ro T . 

CASE Ills The top surface of the enclosure and the mold 
outer surfaces exchange energy by convection 
and radiation with the surroundings, the 
bottom surface being kept adiabatic. 



surface 
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Fig. 4.1. Rectangular enclosure infinite in z-direction. 


4 • 2 . 2 A s s ump t i on s s 


Tlie following assumptions are made in 
the formulation of the problem; 

(i) The physical properties such as density, 
thermal conductivity, specific heat etc, 
o.ro assumed to be independent of temperature 
for the melt, the solid end the mold with 
the exception of the density of the melt 
which contributes to buoyancy forces. 

(ii) The fluid flow developed in the me3_t due 
to thermal gradients is assumed to bo 
laminar . 

(iii) The melt is assured to bo Newtonian. 

( iv) The contact between the solid formed and 
the mold is assumed to be perfect so that 
no contact resistance is introduced. 


4.2.o G- o vo-ming Eg uat ion s ; The equations arc written for 
one half of tlio enclosure because of the symmetry at y= 0. 
The problem may now be stated by expressing the continuity, 


momentum and the energy equations in the melt, and the energy 
equation in the solid and the mold together with the appro- 
priate initial and boundary conditions. Thus in the melt. 
Continuity equation: 
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:.n y-direcTion, 
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ay" p ay 


The first torn in the PJtS of Don. (4.2. 2) is obtained 
by combrr.in •, the vertical pressure gradient and "body force 
term, and a j relying Boussinesq nodol described in Appendix I. 

The pressure terns appearing in Eon (4.2.2) and (4.2.3) 
can be eliminated by differentiating then ~nLih respect to y 
and x, respectively, subtracting and applying the continuity 
eqtiation, Bqn. (4.2.1) to produce the following: equation. 


t 3 u a?) , „ B^u b 2 u .. 3 2 v. „ 3 2 T . 

it ( 'Ty - 1? + "72Ti ; + ' ~'.~J ~ 1 "TIS ~ T 


3Y 


3 X 


3X 3 y 


= £L 


.#1+ v (_lluiL _ 

?.y ay J 


(4.2.4) 

o 

where v is the Laplace operator in x and y coordinates. 
Energjr equation in nelt : 

aT 


3 a 2 ^ 3 T a a 

Tt“ T u "i H v a y 


2 h 2 h 

M— — | + “- 1 — \) (4.2.5) 

ax"" ay 


The problem can be further simplified by writing 
Eqn.( 4.2.4) in terms of vortjcity u and stream function , 
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TThich oases the computational problem 
Thus Lqn. ( 4. 2 .4) becomes 
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(4.2.6) 

(4.2.7) 

(4.2.3) 


The energy equations in the s olio, and the no Id are 
two diner. si ozial conduction equations; 
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(4.2.10) 


These equations, along with the energy balance at 
the interface and the other boundary conditions, are to be 
solved in a domain whose boundaries do not fall along the 
coordinate lines ill any well defined coordinate system. In 
particular, the thickness of the solid formed varies with z 
as we 3.1 as time as the solidification proceeds. This state of 
affairs precludes an analytical solution and requires that, 
even in a numerical solution, appropriations must be made in 
order to obtain results with a realistic computational effort. 
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The initial and boundary conditions for the equations 
discussed above ares 
In the reelt s 

a't = Os T 1= T j _ * (4.2.11a) 

9 I _ 

at y = Ox — y = 0; = 0; y = 0 (4.2.11b) 

at y = l r, s 1^ = T gat ; u = 0* v=C (4.2.11c) 


The 


interface 


ooundary conditions describing the 
can be obtained by writing the heat 


movement of the 
balance equation 


at y=Y* . 


Thus , 

(^(yv) 2 )= pi 4 -H’ 

where 1 = latent heat of solidification. 

The stops involved in obtaiuijy; Eqn.(4.2.11d) are 

discussed in Appendix TV. 

3 T 

at x = Qs —ri = °i a = 0; v = 0 (4. 2. lie) 

d -a. 

at x = is u = 0; v = 0 (4. 2. Ilf) 


Tor cases I and il. 
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For case III, 
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In tlie solid : 
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In the nold: 
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(4.2.12b) 

(4.2.12c) 


(4.2.12d) 


(4.2.12e) 


(4.2.13a) 

(4.2.13b) 


(4.2.13c) 
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(4.2.13e) 
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For case IIx, 
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4*2.4 iT^~diae nsi oncp. izat ion : To make the results more general 
and bo trails form the solution domain into a more tractable shag 
xhe equations are non-dinensi onalized with the help of new 
coordinates x* and y*. The following non-dimensional para- 
meters are defined. 


x*= x/L; y* = y/ Y’ (x,t) ; t* = o^t/D 2 ; 
u*= uD 2 /^!; v*= vD/ a 1 ; *$* = tpD/a^I-; 
o^/a-jL; Tj = (T^T^/Cl.-T^ ; 


T 


* 

s 


f fji _m Wfrp _fp ) • v * = } (lt» t.) „ 

• L r>''^ i i i c ^ y s Y 0 -Y! (z,t7 > 


K = (V^oVdi-V; 4 =^- Y oV(M 0 )i 


(4.2.14) 


The key result of the non-dimensionalizat ion is that 

the range of y* (or y* or y*) extends from zero to one at 
-L s in 

all x* during any time interval. The coordinates x* and y* 
become non-orthogonal as are the transformed coordinates in 
a conventional boundary layer analysis. Further the maximum 
range of temperature between the center line of the enclosure 



and 'tlie mold outer surface "be cones zero to one, during 
solidification. 


% 


Ap pitting these n on-din ensional parameters, Eqs. ( 4.2.5) 

to (4*2.13) get transformed to Eqs (4.2.15) to (4.2.23) (See 
Append in Y); 

Energy equation in melt : 


3 1* a T* 


Yorticity equation; 
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■where Ra = Rayleigh* number 
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and Pr = Prandtl number 
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Stream function equation; 
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Velocity equations; 
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Energy equation in solid? 
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liie initial and boundary conditions become? 
In fixe melt: 


at t 


it 




- 0 : 
= 0: 


T * = 1 ( 4 . 2 . 21 a) 

3T* * 

oj = 0 and v = 0 (4.2.21b) 

3 y x a 7 X 


at y^ = 1: 


T * = l/(l+0); u* = 0; v*= 0 (4.2.21c) 


where 0 = (l i ^ Bat )/(S Bat ^ 0 ) 
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Ste = Stefan lumber 
= I /G (i. - T ) 
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(4. 2,21e) 



at x* = 1 


u* = 0; v*= 0 


Fox cases I and lx. 


(4.2.21f) 


0 


( 4. 2.21g) 


For case III, 


~~~ =~[Bi 2 I* -i- ?.c 2 { (T*+l* ) 4 ~ I'* 4 } ] (4.2.21h) 

B 

where B± 2 = Biot nunber at the top surface of the enclosure 
= ligL/k^ 

he g = radiation constant at the top surface of the 
enclosure, 

= e 2 el 

**0 


In the solid: 


of 7„ = 0: 


at y s = h 


T * = 1/(1 + 0 ) 

O 


T* = I* 


( 4 , 2 . 22 a) 

(4,2. 22b) 


ax x = u: 


(4.2.22c) 


For cases I and II, 


0 


(4.2. 22d) 
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For case III, 


Li§' = -(^)[3i/ s + So 2 {(V<) 4 < 4 ,] 
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In the mold; 


(4.2. 22e) 
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(4.2.231) 


ior case 


1 * = 0 


(4.2.23c) 


For cases II and III, 


[Bi,T* -i- Re, { (1* + 1* ) 4 -l* 4 > ] 
1 m l v m c c /J 


(4.2.23d) 


trhere Bin = Biot number at the nolcl outer surface 




Lc^ = Radiation constant at the mold outer surface. 
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Fox cases I and II, 

3 1 * 
m 


For case III, 




(4.2.251) 


(4.2.25s) 


4.5 SOLUITOh METHODOLOGY; 

-Hie crus of the problem lies in solving the simul- 
taneous, partial differential Eqs (4.2.15) - (4.2.20) for 
the boundary conditions given in Eqs (4.2.21) - (4.2.25). 

It is evident that the problem cannot be solved analyt i c ally 
because of the complexity and the non-linear nature of the 
equations. This leads to the obvious choice of numerical 
methods. 

In the computational model it has been assumed that 
there is a small time lag bet ire en the heat withdrawn at the 
interface and the resulting interfo.ee motion. Specifically, 
during- a smell interval of time, the fluid floir in the melt 
and the heat transfer in the noli, solid and the mold i-rere 
solved assuming interface to be stationary. The heat trans- 
ferred to the interface during thin interval was then used to 
compute the snail finite displacement associated with the 
solidification. With the new domains for the melt and the 
solid, which correspond to the displaced interface, the 
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solution for the next time interval was performed as just 
described. The nature of Eqs (4.2.15) - (4.2,23) implies 
that the methodology of solving them can be divided into 
four sections, viz, 

(i) Solution of vorticity and energy equations 

which are parabolic in nature - .^ch aspect to time * 

(ii) Solution of elliptic stream function 
equation. 

(iii) Boundary conditions for the vorticity 
equations. 

(iv) Solution of moving interface equation. 


4.3.1 Vor ticity and Energy Equations i The. p arab oli c , v or t i~ 
city and energy equations were solved using the Alternating 
Direction Implicit (ADI) technique introduced by Peaceman and 
Jiachford (1955). This method was particularly chosen as it 
has been found suitable for s olving similar kinds of problems 
by Wilkes and Churchill (19 66) and, Szekely and Todd (1971). 
Also known as the method of variable direction, the method 
makes use of splitting the time step to obtain an implicit- 
method which requires only tbe inversion of a tridiagonal 
matrix. The non-linear convective terms were written in 
their central difference forms as both the forward and 
backward forms led to numerical instability. Por high j-'-ayleigh 
numbers (Ra >_ 10^) the convective terms were modified 

using 'Second upwind differencing method' (see Appendix VI) 



as the program became numerically unstable# The tridiagonal 
matrix formed was solved using ‘General iridiagonal Algoritlrn’ 
discussed by P.oaclie (1976). The accuracy of the grid size and 
the tine stops was checked by doing the calculations for 
different grid sizes and tine step . For example, computations 

rz 

were done for ha = 10, by choosing grid sizes of 10x10, 12x12, 
15x15 and 20x20, The variation in the results for the first 
three cases was less than 1 per cent. In the last case, hie 
variation is about 5 percent which is due to round-off errors. 
The computational tine depends upon the grid size chosen. 

For ila = 10^, when the grid size was 10x10, the CPU tine 
was 180 sec. for t*= 0,27; for Fa = 10 6 , the grid size was 
50x50 and the CPU tine was 900 sec. upto t*= 0.18. An example 
of the All algorithm written for the velocity Eqn. (4,2.16) 
is shown in Appendix VII. 

4,3.2 Stream Fun ctio n Equation ; The elliptic stream function 
Eqn* (4.2.17) was solved using Successive Over-Relaxation 
(SOR) method discussed by Forsyth and Wasow (1969), This 
method was particularly chosen as the convergence rate was 
faster than the other methods such as Southwell's residual . 
relaxation method. Gauss— Seidel iteration method euc. 

Eqn. (4*2.17), when written in finite - difference foir 


gives. 
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(4.3.1) 


This can he simplified to give, 
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(4.3.2) 


where 

n=<vi> 2 1& 2 

y 2 = ( f ,) 2 

i,j denote space while k denotes iteration step. In hq_s. 
(4.3,1) and (4.3,2) the superscript * is omitted in 4 > , w , 

Ax mad A y for convenience. Convergence is assumed when 
, k+1 • k+1 . • 

I * . . - . I = 0.001. 

i,3 1,3 
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^ 0 p^. is the optimum relaxation fortor, the value of 
which is to he found by trial and error, for this parti- 
cular problem ft . = 1 . 5 . 

^ opt 

4.3.3 Yorticities at B oundaries : 


Special attention has to he given to the vorticities 

at the boundaries as these values are not known explicitely. 

If P . is the stream function at the wall, 

w, J 

expansion of 4 . by the Tailor series g:h r es, 

W , J 


^w- 1, 3 
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AX + ~ 
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+ 0( A x)' 


(4.3.3) 


and. 


3 4 

TF 


j = v | =0 (no slip condition) (4.3.4) 

w,j w,j 


4 . = 0 

w, 3 


Thus Eqn.(4.3.3) becomes, 

1 3 2 4 ( . v 
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(4.3.5) 


(4.3.6) 


from the stream function Bqn.(4.2,17) 


( D ^2 34 | 

w. T s = - l r ) 2 ' 
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Subst ituiting Eqn.(4.5.6) in Ban. (4.3.7), 

2 4. 

/ U \ ^ 
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( 2) 2 


W-1,.1 


(4.3.8) 


( A x)‘ 


^w-1,3 "^ ie i n ^ er i° r stream function discussed 
in Section 4.3.2. The vorticity values calculated using 
higher order terns of the Tailor series were found to create 
stability problems. This is in confirmation with Samuels 
and Churchill (1967) who have reported similar problems. 

4.3.4 Moving Boundary Ecmation ; The moving boundary Eqn. 
(4.2.21d) is a cubic polynomial similar to the one described 
in Chapter III. This has been solved using Chebyshev 
approximation discussed by Carnahan, luther and Wilkes (1969). 

4*3.5 Choice of G-rid System ; The choice of grid system is 
very important in two dimensional solidification problems 
as the interface is moving with respect to time, the shape 
of which is not known. In most of the conduction -based 
problems Cartesian grid systemhas been chosen for computation. 
Figure 4.2a illustrates Cartesian grid of points formed 
between the interface and the other boundaries in the melt and 
the solid. This may be regarded as the most obvious choice 
of the grids; yet it has the following disadvantages: 

(i) It puts only a few points in the melt wherever 

that is thin (i.e. at the bottom of the enclosure) 
and so c ann ot provide high accuracy there. If more 
points are provided in those regions, it increases 
storage space and compu.tation tine enormously. 

( ii) Calculation of the grid length near the interface 
is complex. 



Interface Irterf ice 



,4-2a Cartesian arid system Fig.4-2b Grid system in the 

present work 



(iii) Very few grid points con be provided in the 

solid region during initial stages of solidi- 
fication. 

This led to the choice of a new grid systen described 
by Patankar and Spalding (1567). Figure 4.2b illustrates the 
grid system adopted for the present problen. In the a* coordi- 
nate, the lengthof the enclosure has been divided into grids 
of equal size. At any instant of tine and for any x , tne 
width of the licit or solid region has been divided into equal 
number of grids. This method of division of grids is devoid 
of all the disadvantages mentioned for the Cartesian grids and 
has been successfully used by Spalding (1977) in solving boundary 
layer problems. Spatial grids of 10x10 were established, each, 
for the melt, solid and mold. For Bn. ^10^ closer grid was provi- 
ded to take care of the stability problems. 

4.5,6 C or in ut at ion a 1 Stops ; The following stops were adopted in 
solving Bqs.(4. 2.15) - (4.2.23). 

****** 

(i) At tine t =0, the initial values of T-^, ^s’ m’ W , * 

and v* were established. 

(ii) At any intermediate time, thickness of 'die solid iorned 
was calculated through Eqn.(4,2. 2ld) . The temperature gradients 
involved in Eqn. (4.2.21d) were calculated using the tenperaruro 
fields known from previous tine steps. The field variables 

■** w * u* and v* were redefined for every time stop to take 
cognizance of the change in the relationship between the physical 
and transformed plane, At any time step, the re -evaluation ox 



these variables in the grid points was done by interpolating 
through the j£nown values of the previous tine step grid points 
using Lagrange formula of interpolation. To conserve the 
momentum. the interpolation of tbe variables was perfornod by 
first evaluating the stream. function at the locations of the 
interest and then differentiating it. Sinilarly the tenperature 
was interpolated through the evaluation of the enthalpy. 

(iii) For the known values of velocities, the tenperature was 
advanced by one tine step by solving Eqn.(4.2.15) using ADI 
technique . 

(iv) From the tenperature, estimated through step (iii), inter*- 
ior vorticity field was calculated through Eqn.(4.2.16) using 
the ADI technique. The boundary vorticity values conputed for 
the previous tine step were used for the calculation. 

(v) The stream function Eqn.(4«2. 17) was thaa solved using 
SOB. method, 

(vi) The boundary vorticities were then calculated using the 
stream function conputed via step (v). 

(vii) Eron the stream function, the velocity components were 
calculated. 

(viii) The last step in the computation involved the calculus i.~. 
of tenperature fields in the solid and the mold through Bqs. 
(4.2,19) and (4.2.20), respectively. 

The computational cycle was repeated until desired 
time or desired amount of solidification was reached. 

Schematically, the computational procedure is illus- 
trated in Figure 4.3, A listing of the computer program used 
on DEC 1090 computer is available with author. 




Start 


Read Data 


Initialization of temperatures, vorticity 
stream fu net i on . velocities 


Calculate new t 
Solid t h ickness ( 1 -6 ) 


Temperatures in melt reg ion , eqn.(4.2.15) 


Vorticities at interior 
points, eqn. (4-2-1 G ) 


Iterate for new j/r 
at all points, eqn. (4.2.17) 


Vorticities at the 
Boundaries, eqn. (4. 3. 8 ) 


Velocit ies u . v 
eqn, (4.2-18) 


Desired ^ 
/time reach' d 
or 

desired sol. id if teat ion 
; d c h i Vfed • 


Temperatures in solid 
and mold , eqn (4-2-19M4-2-20) 


chart i 




4.4 RESULTS AH I) DISCUSSION; 


The following results are presented; 

(i) Isotherm and streamline patterns developed in 
a rectangular enclosure, in which the melt does 
not solidify at any time. (Figures 4.4 - 4.6). 

(ii) Temperature and velocity patterns and interface 
movement developed in a rectangular enclosure in 
which solidification occurs. The top and bottom 
of the enclosures are kept adiabatic and mold outer 
surface is kept at a constant temperature. 

(Figures 4.7 - 4.16). 

( iii) The interface movement during solidification 

when the top and bottom surfaces are adiabatic 
and mold surface exchanges energy by convection 
and radiation. (Figures 4.17 ~ 4.19). 

(iv) Isotherms, streamlines and the interface patterns 
developed during solidification when the top 
surface and the mold outer surfaces exchange 
energy with the surroundings by convection and 
radiation, the bottom surface being kept adiabatic. 
(Figures 4.20 - 4*28). 


4.4.1 Thermal and F luid Flow, jin M el^^tho ut . ^ojm^ ^icahion 
This case is considered mainly to have a chock, on the 


comp ut o.t i ona 1 program wr i b u cr- 
on d to get an idea above; the 


for the solidification case 
flow pattern that would be 



developed due to temperature gradients 


Figure 4.4 represents the transient streamlines 
and isotherms when the top and bottom surfaces of the 
enclosure arc kept adiabatic while the sides arc- kept at a 
constant temperature. Figure 4.5 shows the streamlines and 
isotherms when the top surface is cooled by convection and 
radiation, sides at a constant temperature, the bottom 
being adiabatic. The hot fluid rises up near the center of 


the enclosure and the cold fluid comes down near the side 
walls whore the temperature is minimum. The flow pattern 
for Fa = 105 reveals a rolling vortex centered at approxi- 
mately x*= 0.5, y*= 0.5. The maximum value of the stream 
function 4 , occurs at the center of this vortex. As tne 
Rayleigh number increases the vortex center moves to the 
right -hand side. At very high Rayleigh number (Ra = 10 ) , 
the entire enclosure is more or less occupied by the ris.ng 
hot fluid. The streamlines near the wall move closer to- 
gether as a result of high flow in this region. 

At low Rayleigh numbers, the isotherms are almost 
vortical indicating the predominance of conduction. Tno 
isotherms take on inward turn at the top when the enclosure 

is cooled from the top, as shown in Figv.ro 4.5. As the 

n rt/vf'h rv^-mc* ripvTf’tc 'from tlie vertical 
Rayleigh number increases, lsotflOiim, novice n»* 

position as convection starts to play a major role. 

Figure 4.6 represents the streamlines pattern when a 
fluid of Pr=l, is heated in a rcctau C ular enclosure. A 




Fig.4-4 Streamlines and isotherms for different Rayleigh numbers at . 
t° s 0*027- without solidification, top and bottom surfaces adiabatic 
sides at a constant temperature 
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m >i«i mi 


Ra = 5070 
Pr =1-0 
t* = 0-05 


At \%0,T*--Q 

At t >0,17 = 0*5 at y =Y 0 ( he ati ng ) 

Present work 

• Szekely and Todd (1971) 



0 025 T v ~0 0-- c To 

y/Y 0 

Fig. 4-6 Comparison of the streamlines of present 
work with liter aturc-wi t hout solidification 
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comparison of tho results, obtained from the present 
numerical scheme -with, the results reported by Szekely and 
Todd (1971) shows a reasonably good agreement . 

4.4.2 The rmal and Fluid Plo w . in lie It -with Solidification , 

Top and Bottom of the Enclos u re Ad i abatic. Hold Wall 
at a Constant Temperature ; The temperature profiles 

in the molt, the solid and the mold at various heights, with 

non— dimensional tine t* as the parameter are shown in 

Figure 4.7. As the tine progresses the temperature decoys 

in the nolt . Tho temperature of tho molt is found to 

increase with the height because of the natural convection 

developed in tho melt. Initially, tho temperature variation 

with the height of the enclosure is very little since 

conduction is dominating . Subsequently, natural convection 

currents are developed due to thermal gradi ent s , which in 

turn take the hotter melt to tho top of the enclosure driving 

the cooler melt to the bottom. 

The time dependent velocity profiles in the melt are 

shown in Figure 4.8. The velocity profiles arc symmetric 

about y*= 0 and the net upward flow equals the not downward 

flow. With the progress in solidification, the vortical 

velocity of tho nolt increases, roaches a maximum value at 

some intermediate time and then decays. At time t*= 0.03, 

tho convection currents have just set in. The velocity 

near tho top (x/L=0.9) and near tho bottom (x/L=0.l) of the 

enclosure is found to be the same, and the melt-solid 



Solidification 

temperature 


Ra . 

5x 10 3 

Pr ^ 

0-98 

Ste =■ 

0-5 

0 -- 

1-0 

L/Y 0 = 

2-2 



0-2 0-4 0-6 0-8 1 - 

y /Y o 

ig.4.7 Temperature distributions in the three regions 
at different X ^ f or constant temperature at 
the outer surface of the mold and adiabatic 
at top and bottom of the- enclosure. 




Vertical velocity 



y/Y 0 

Fig. 4. 8 Vertical velocity distribution in the melt 
at di tferent X* 
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interface is almost parallel to the wall. As solidification 
progresses, the slope of the interface changes, and the 
velocity near the top becomes more than that near the 
bottom surface of the enclosure. A comparison of Figures 4*7 
and 4.8 at ecjual heights but different times indicates that 
the thermal field and the fluid motion decay almost at 
equal rates. This is true since Pr=0.98, which implies 
that the thermal and momentum diffusive. t ios are nearly 
equal. 

The streamlines and isotherms for the Rayleigh 

W Z 

number 5 x lCr and at time t*= 0.27 are shown in Figures 4.9 and 
and 4.10 respectively. The flow pattern during early stages 
of solidification is very nearly similar to the one reported 
in Section 4.4.1 and is not reported here for brevity. From 
Figure 4.9 it can be easily observed that the general 
characteristic of the flow pattern is upflow near the center 
of the enclosure, where the melt is hotter and downflow near 
the melt -solid interface where the molt is cooler. Both tho 
streamlines and isotherms indicate that at larger tine, the 
thickness of the molt layer is greatest at tho top and least 
at the bottom of the enclosure. It is interesting to note 
that the streamlines of higher values arc restricted to the 
top portion of the enclosure where the temperature and 
hence the flow is maximum. At high Rayleigh number (Ra = 

5 x 10^) secondary cells are observed near the center of 



if/ - if i ^/pnax 


Pr - 0-98 

St<2_ 0-5 
0 --V0 
L/Y 0 = 2-2 

k s /k m = 1-2 


0-25 


0-50 


0-75 


Fig. 4.9 Streamline pattern in the melt forRa=5x10 
and t*=0-27, for constant temperature at 
the mold surf ace, adiabatic at top and 
bottom of the enclosure. 
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the enclosure during initial period of solidification as 
shown in Figure 4.11. 

An examination of Figure 4.10 would indicate the 
predominance of the natural convection there the isotherms 
are non-linear, contrary to the linear vertical isothorns 
one would expect in the co.sc of pure conduction. Isotherms 
of higher values are concentrated near the top surface of 
the enclosure where the temperature is maximum* 

The dependence of the interface movement on the 
Rayleigh number is shorn in Figure 4.12. At low Rayleigh 
number (Ra = 5 x 10 ) where the natural convection is minimum, 
the shape of the interface is found to be almost vortical 
and parallel to the mold. As Ra approaches zero, the inter- 
face will be exactly vertical since conduction is the only 
node of heat transfer then. As Ra increases heat transfer 
duo to convection increases and the slope of the interface 
changes (Figure 4.12). At low and moderate Rayleigh numbers 
the computation is carried upto non-dimensional tine t = 0.30, 
whereas for high Rayleigh numbers, the computation is res- 
tricted to a time t*= O.tg, due to enormous computation 
time required in the latter case. 

The movement of the interface with tine if for 
different aspect ratios ( L/ 2Y =1.1 and 5.5) is shown in 
Figure 4.13. At low aspect ratio, where the convection 
is maximum, the slope of the interface is high and it increases 



Fig. 4-10 1 sotherm_pattern in the melt for 

Ra-5x10 3 and t*=0'27,for constant 
temperature at mold outer surface 
adiabatic at top and bottom of 
the enclosure 




Pr -0-98 ; <p = 1-0 
Ste-0-50 ; L/Y 0 -_2-2 
- 1-2 
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0-50 

y/ v 0 
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Fia 4-11 Streamline pattern in the melt for 
y Ra r 5 x 1 0 5 and t*=0-03,for constant 

temperature at mold surface , adiabatic 
at the top and bottom of the enclosure 



= 5x10 
5x10 

5x10 


t* 0 


Pr -0-98 

St <2 -0- 5 
0 -10 

L/Y 0 _2- 2 
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Fig. 4-12 Dependence of the interface movement 
on Rayleigh number at different time, 
for constant temperature at mold 
surface , ad iabatic at top and bottom 
of the enclosure 




Ra = 5xl0 3 ; Pr zO-98 
St<2 - 0- 5 ; 0^1-0 
k s /k m -.l.2 


ig.4-13 Interface movement at various time intervals dor 
L/D ratios of 2-0 and 10*0, for constant temper an 
at the mold surface/adiabatic at the top and 
bottom of the enclosure 
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with. tine. At high aspect ratio, the interface is almost 
vertical except near the top and "bottom of the enclosure, 
where the boundary effects start influencing the shape of 
the interface. Ramachandran, Gupta end Jaluria (1981b) 
suggest that one can assure- the heat transfer to "be uni— 
directional when the aspect ratio is high. 

The influence of Prandtl number (Pr) on the interface 
movement is shown in Figure 4.14, As Pr increases the 
interface slope increases ujito Pr=10. Put erostinglj , Prandtl 
number greater than 10 does not have any influence on the 
interface movement and an asymptotic solution is reached. 

This supports the results obtained by Sparrow, Patanhar and 
Ramadhyani (1977). This beh.aviour can be better understood 
from the vorticity Eqn.(4. 1.16) , which can bo simplified for 
high Prandtl number as: 


(B) 2 + (|h~ l ~3 = (xKfrH-?-) 3 2a — • i (4.4.1) 

n a x*^ 1 ayj 1 x o syj 

It is important to note that the solution of Eqn. (4.4.1) is 
much simpler than the parabolic vorticity Eqn. (4. 2. 16 ) . 

The effect of St of an number (Ste) on the rate of 
interface movement is shown in Figure 4.15. The melt-solid 
interface is found to proceed faster when Ste is small, as 


expected. 

The dependence of the solidification on superheat (0) 
is given in Figure 4. 16, where the thickness of the solid 



Fig. 4-14 Dependence of interface movement on 
Prandtl number at different times, for 
constant temperature at the mold surface, 
adiabatic at top and bottom of the 
enclosure 




y/Y 0 

Fig. 4*15 Dependence of interface movement on 

Stefan number at various time intervals 
for constant temperature at the mold 
surface, adiabatic at top and bottom 
of the enclosure 





Fig. 4-16 The effect of superheat on solid formed, for constant 
temperature at the mold c urfuoe, adiabatic at top 

and bottom of the enclosure 




formed is plotted against time with, f) as the parameter. 

The rate of movement of solidification front decreases with 
increase in 0. This is expected from the results of one 
dimensional analysis, as large J# implies higher value of 
superheat and hence a decrease iri the solidification rate. 

4.4.3 Thermal and Fluid glow in H olt -with Solidification ,. .Top 
and Bottom Surface A di abatic, and Convection. and 
Rad ia tion at the Hold Outer Surface ; The variation 

of temperature and the vertical velocit; with height at the 

center of the enclosure and the dependence of interface 

movement on Biot number (Bi^) and Radiation constant (Hc^) 

ar o discussed in this section. The streamlines end isobliej.m 


patterns and the dependence of interface movement/ on non- 
dimensional parameters such as Ra, Pr, Ste and aspect ratio 
differ only quantitatively from the previous case (Section 
4.4.2), olid are not reported here. 


The variation of temperature with height, at the 
center of the enclosure and at the mold outer surface are 
shown in Figure 4.17. Initially (t*= 0.06) there is very 
little variation in the center line temperature in the in- 
direction, as conduction is dominating. At larger times 
the buoyancy forces play a predominant role and the temperature 
increases with height. It is also found from Figure 4.17 
that the mold outer surface temperature does not vary with 
the height at any time. This is true since the mold thickness 
is very small compared to the height (l/d=20). This suggests 
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ig.4-17 Temperature profiles at the Center lint of 
the enclosure arj mold outer surfaces at 
different times -Top and bottom of the 
enclosure a ii abatic , mol I surface 
converting an] raliating 



that one car, safely assume the heat transfer in the mold 
to be •unidirectional if the L/d ratio is very high. 

Iho variation of the center line vertical velocity 
u with the height is shorn in Figure 4.18. The curves are 
symmetrical about x/L=0,5 and the nanism velocity occurs 
at x/!=0.5. 

Figure 4.19 shows the dependence of interface move— 
moni/ on Biot number and Radiation constant. An increase in 
Biot numb er or Radiation constant increases the rate of 
solidification. It is seen from Figure 4.19 that very high 
Biot number (Bi-^ > 15) or Radiation constant (Rc^ > 500) does 
not have any influence on the rate of solidification and an 
asymptotic solution is reached. This is true as high Biot 
number or Radiation constant indicates a high convective 
or radiative cooling and the mold approaches the ambient 
tonpe nature. 

4.4.4 T hermal and Fluid Flow in .Melt-with Solidification , 

Bottom Sur face Adiabat ic, and. Convection and Radiation 
at the Ton Surface and the Hold Outer Surface : The 

computed results are presented mainly to study the effects 

of non-dimensional parameters on the rate of solidification, 

streamlines and on the isotherms. Of the non-dimensional 

parameters, attention is nain3.y focussed on the Rayleigh 

number, the Biot number and the Radiation constant. Rayleigh 

number is particularly chosen as it characterizes the role 

of natural convection. Ramachandr an, Gupta and Jaluria( 1981c) 
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Fig. 4*18 Velocity profiles at the c^nLr line of the 
enclosure for different ti men -Top and 
bottom surfaces a liul itic, mol i surface 
convccting an 1 radiating 



Fig. 4-19 Dependence of interface movement on Biot 
number and Radiation constant-Top and 

bottom surface adiabatic, mold surface 
convectinq and radiating 
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have discussed the importance of Rayleigh number in 
characterizing the role of natural convection during two 
dimensional solidification. 

The influence of the Rayleigh number on the interface 
movement is shown in Figure 4.20, Initially, for low 
Rayleigh number the shape of the interface is vertical and 
parallel to the mold except near the top of the enclosure, 
where it flares in due to cooling at the top surface. As 
Ra increases the interface deviates from the vertical 
position, slopes outward and then flexes inside duo to local 
cooling from the top. This effect is quite marked at tine 
t*=0.18, when the natural convection has started playing a 
significant role. This is contrary to the pure conduction 
where the interface slopes inside at any point inside the 
enclosure. At time t*=0,36, the effect of cooling at the 
top surface is quite significant which changes the slope- of 
the interface gradually inside even at the middle of the 
enclosure. For high Rayleigh number (Ra = 1x10 ) computer ion 
is restricted to t*= 0.18 due to enormous computational time 
required because of small time srep. 

The role of the Rayleigh number can bo bettor under- 
stood by computing the heat transferred across the interface 
for different Rayleigh numbers and comparing with the 
results of pure conduction case. The quantity of heat 
transferred across the interface of length ax and unit wio-bh 
in the z— direction can bo written as 



x 10 
x 10 
x 10 


L / Yo = 2 • 2 


Fig. 4-20 Dependence of interface movement on 
Rayleigh number at different times - 
Bottom surface ad iabatic , top surface 
and mold outer surfaces convecting 
and radiating 
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On non-dimonsionalization this becomes, 
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rlierc , Q = °C LY (f.-T ) 
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Equation (4.4.4) is integrated for every time stop 

and tho value of Q/ Q q is found for different Rayleigh numbers 

and for the pure conduction case. Tho ratio Q, ,, \/Q 

- { cond+ c onv ) c ond 

is plotted against tine as shown in Rig-are 4.21, for different 


P 


is the heat transferred at the interface 


lay 1 c igh ironb e r s . 

^( cond+conv) 

with consideration given to the natural convection in the 

nelt apart from the heat transferred by contfaction. 

Q , is tho heat transferred bv conduction at the 

cond 

interface without considering the natural convection in the 
melt. As solidification progresses tho ratio Q ( C0nd+C0OT )/Q 001l4 
increases, attains a naxinuri and then decreases. In case of 
low Raylc-igh numbers the occurrence of these events is shifted 
to larger times which is not shown in the figure* Initially 
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the ratio increasos as appreciable heat is transferred by 
natural convection, which is quite predominant in the case 
of high Rayleigh numbers. As solidification progresses the 
thickness of the melt region decreases and so does the 
heat transfer by natural convection. This explains the 
decrease in heat transfer after some time. 

Figures 4,22 and 4.23 slow the dependence of tie 
interface movement on Radiation constant and Biot number 
r e s p e c t i vc ly , It is found that an increase in Radiation 
constant (or Biot number) increases the rate of interface 
movement, which is evident as discussed in Section 4.4.3, 
for different set of boundary conditions. 

The streamlines for the Rayleigh number 1x10^ at 
t*= 0,06 and t*= 0,30 ar^ shown in Figures 4.24 Mid 4,25 
respectively. Figure 4.26 shows the streamlines for the 

jr 

Rayleigh number of lxl0 ? at time t*= 0.06. The streamlines 
are normalized with respect to the maximum value of the 
stream function which is shown as a dot in all the figures. 
The interface, center and the top of the enclosure represent 
the streamline 4* =0, The direction of the motion of the 
melt along the streamlines is denoted with small arrowheads. 

For Ra=10^ at t*= 0.06, the streamlines arc almost ^ 
symmetrical about x/L=0#5. This is because the interface 
is vertical and parallel to the mold wall except near the 
top of the enclosure. In this case the flow in the raelt is 
similar to the flow that would develop in a rectangular 









Dependence of interface movement on 
Biot number at different times - Top 
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cavity in which, no solidification occurs. The streamlines 
move closer together near the interface as a result of 
higher flow in this region. As solidification progresses 
the shape of the interface deviates from the vertical 
position and the streamlines are no more symmetrical about 
x/ L=0.5 . file melt moves upward near the center ox the 
enclosure, flows gradually towards the interface as the 
temperature decreases with height near the top of the 
enclosure. The melt then flows down against the retarding 
effects of buoyancy. As tho molt is slowed down it turns 
inward and risc-s to a height determined by its buoyancy , 

A comparison of Figures 4,24 and 4.26 indicates, for the sano 
time as Rayleigh number increases the no It moves towards the 
interface even near the center of the enclosure. She naximun 


valuo of tho stream function is increased by more than ten 
times and this value occurs very near to tho bottom of tho 
on clo sure . 

Figures 4.27 and 4.28 show the isotherms for time 

3 5 

t*= 0.06 for Rayleigh numbers 1 s 10 r and 1x10 respective^. 
In the case of R.a=10^ where the flow is yet to develop the 

5 

isotherms almost taiio the shape of tho interface. For Ra=10 
the magnitude of tho flow is quite high and tho shape of tho 
isotherms of higher valxtcs deviates completely from the 
interface shape. Ac the top surf are is cooled by radiation 
and convection there is a 


local cooling of the nc-lt near the 
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Fig. 4-27 Isotherm pattern in melt for Ra = 1x10 
and t*=0-06- Bottom surface adiabatic, 
top surface and mold outer surfaces 
con vec ting and radiating 
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Fig. 4-28 Isotherm pattern in melt for Ra =1x10 
and t*= 0*06 - Bottom surface adiabatic 
top surface and mold outer surfaces 
convectinq and radiating 
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top of the enclosure and the warmer portion gets stratified 
in the center in between the cooler portion of the nelt. 

4.5 COHCLUSIOITs 

An analysis lias boon carried out for the s olidifi cation 
in a rectangular enclosure. The equations are written for 
the temperature, vorticity, stream function and velocities 
in the nelt. In the analysis, various combinations of 
boundary conditions are considered at the top and bottom 
surfaces of the enclosure and at the mold outer surface. 
Consideration has also been given to the heat conducted 
through and stored in the solid foruod and the mold. The 
coupled, non- linear, simultaneous equations arc n 011-din on si du- 
alized and solved using ADI and SOU finit c-diff orencc 
techniques. The thickness of the solid formed at every 
time step ha.s been found by solving the energjr balance 
equation at the melt-solid interface. 

The temperature and velocity profiles are plotted for 
various heights in the enclosure. The temperature is found 
to increase with height when the top surface of the enclosure 
is kept adiabatic. In the case of radiative and convective 
cooling from the top, the temperature first increases and then 
decreases with height. It has also been found that the 
natural convection has significant effect on the shape of the 
interface. The vertical velocity component and the temperature 
decay at almost equal rates when the Prandtl number is close 
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to ■unity. The effect of natural convection on solidification 
is cleanly shown in the streamlines and isotherms plotted 
for different Rayleigh numbers and for different "boundary 
conditions. Finally, the rate of solidification is also 
studied for various non-dimensional parameters such as 
Rayleigh number, Prandtl number, aspect ratio, Stefan number, 
Biot number, Radiation constant and the extent of superheat. 



OHA gggt 5 


EXpaiiTiuHTs 

5.1 lifi'RODUCTIOI 

The experimental study of solidification of melt 
in an enclosure of finite dimensions lias been made by many 
investigators. Host of the experiments have been conducted 
to study the shape of the interface and the temperature 
profile in the melt and the solid wnthout considering the 
natural convection in the melt. The main reason for this is 
the sensitivity of the experiments to different parameters 
involved, as the results can vary to a large extent from one 
experiment to another depending upon the extent of superheat 
in the melt, the exact position of the mold, the cooling 
rate at the mold surface and the attitude of the enclosure with 
respect to gravity. In fact only recently attempts are 
being made to compare the experimental results of solidifi- 
cation until the analytical work. Sparrow, Ramsey and Kemink 
(1979) have conducted experiments in an annulus of concentric 
cylinders and compared the results with the theory developed 
by Sparrow, Patankar and Ramadhyani (1977). 

In most of the experiments, a low melting paraffin 
or ammonium chloride solution has been used as the solidi- 
fying material. A low melting paraffin is mainly used 
because it freezes at a constant temperature. Ammonium 
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chloride solution is very useful to study the solidification 
of binary dendritic systems because dire to its low entropy 
of fusion the material solidifies very much like metallic 
alloys. In most of the cases only temperature measurements 
using thermocouples have been made as it is very difficult 
to measure the fluid flow developed in the melt. 

In this chapter, the experiments performed for the 
solidification in an enclosed region with paraffin wax 
(properties in Appendix VIII) as the solidifying material 
are described. Paraffin wax (melting point 60.5°C(i40.9°P) 
is particularly chosen because of its low melting, easy 
handling nature . Water is used as a coolant for solidifi- 
cation. The temperature of the solidifying wax is measured 
with thermocouples kept at different predetermined levels 
inside the cm closure. The recording of the temperatures 


has been done with the help ot a programmable Data Acquisition 
System. The variation of temperature with tine is plotted 
and compared with the theory. Prom the temperature vs time 
plots, the time taken for the melt -so lid interface to touch 
a particular thermocouple is interpolated. The position of 
the interface is then plotted against tine and compared with 
the analj-tical solutions given in the previous chapters 
(Chapters 3 and 4). The experiments are conducted for 
different aspect ratios and superheats of the melt, in the 
laminar range of natural convection. 
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5.2 EXPEItll-SmAL APPAItATUB s 

The heart of the experimental sot up is the tost 
cell where tho solidification takes place. In addition to 
the test coll, tho apparatus includes auxiliary systems 
for controlling tho temp eraturo of cooling rater* Data 
Acquisition System to record temperatures and a heating 
device for melting the xr<xx. 

Separate tost colls were fabricated for conducting 
exp er incut s of one dimensional and two dimensional solidifi- 
cation cases. Figure 5.1 shows a schematic diagram of the 


experiment al set up used. 

5.2.1 lost Coll for one Dime nsional Case ; The test cell 
for tho one dimensional experiment consisted a plexiglas 
rectangulo.r container of dimensions 14r50x50 cm ( 5 .5x19 .7x19 .7 
in), into which two 0.64 cm (1/4 in) thick aluminium plates 
were inserted, thus dividing it into three compartments. 

Tho paraffin wax was solidified in tho central compartment 
of dimensions 3.0x50x50 cm (5.9x19.7x19.7 in) and the cooling 
water was circulated in the side compartments of dimensions 
2x50x50 cm (0.79x19.7x19.7 in). T'lic aluminium plates formed 
tho two vertical opposing walls of the mold, A Plexiglas 
sheet of 1.25 cm (1/2 in) thick, used to cover tho molten 
wax in the centra], compartment, formed the top surface 
of the enclosure. The top surface was kept at a height of 
47 cm (18.5 in) from the bottom of the enclosure (l/2Y 0 =4.7). 



Data 
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Ail opening was provided, at the top surface for pouring 
the molten wax into the central coup art: icnt . The wax can 
he drained through an outlet kept at the bottom of the 
central compartment. All the- sides, top end the bottom of 
the enclosure were insulated with 1.25 cm (1/2 in) thick 
asbestos and 2.5cm (1 in) thick styrofoam shoots. 


5.2.2 lest P ell .for t wo Dimensional Case ; The test cell 
for the two dimensional experiments was similar to the one 
described above with slight modifications. The sis*c 'of the 
cell was much smaller than the previous case so that the 

A 

Rayleigh number was well within the- laminar range (Ra^SxlO ). 
The test cell., with dimensions 2x20x2ocm (0.79x7.9x7.9 in) had 
two walls of stool (20 x 20 cm; 7. 9x7* S' in ) and two side 
walls of Plexiglas (2 x 20 cm; 0.79x7.9 in ). The steel 
walls were 0.2 cm (0.079 in ) thick and acted as mold as 
considered in the theoretical analysis (Chapter 4). The 
stpol walls had jackets (5x20x20 cm; 1.97x7.9x7.9 in ) in 
the ' other sides for circulating cold water through them at 
controlled t emperatur o . 

The height L to which the molten wax was poured 
can be varied by placing 1.0 cm (0,4 in ) thick Bakelito 
sheet in the test coll, at the- levels of 5/5, 5,6 cm (1.38, 
2.0, t 2.4 in ) from the bottom of the container. Grooves 
were made a„t these levels for sliding the Bakelito sheet 
which served as the top surface of the enclosure. 'The molten 
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wax can. Vi drained out through, an outlet provided at the 
bottom of the test coll. 2 ho whole tost coll was insulated 
with 2.5 cm (1.0 in ) thick styrofoam sheets. 

5.2.5 T emp or atur o M o assure m on t s : Temperature of the soli- 
difying wax was measured using a double insulated, 24 gauge, 
Iron-Constant an thermocouples. The thermocouples wore 
sheathed inside copper capillaries, fixed inside the tost 
cell. Eighteen thermocouples 'were used for measuring the 
temperature of the wax, kept at six different levels in the 
^-direction, each level containing three therno couples. 

Separate thermocouples were provided for measuring the 
temperatures of the mold wall and cooling water. For one 
dimensional experiments, only eight tlx or ao couples were used, 
all kept at a height of 25 cm (9.8 in) from the bottom of 
the test coll. Figure 5 . 2?b shows the photograph of the 
thermocouple assembly used for the two dimensional experiments. 

The output of the thermocouples was recorded by a 
5052A Hewlett Packard, Data Acquisition System (DAS). The 
DAS is basically a high speed automated instrument that 
provides an on-line analysis for the real time experiments. 

The system is composed of 3 45 5 A high resolution, digital 
voltmeter (Dili), 3495A scanner, 59309A real tine clock and 
9825A programmable calculator, all interfaced through 98034A 
HP -IB interface. 






Tlio thermocouples iroro connected to different 


channels of the scanner which is ba.sica.ll3r a, pr ogramnab lo 
relay bon. l'lirotigh the scanner relays, therno couple output 
signals were routed to JNll at specific time as controlled 
by the calculator program. She speed ox the scanner is 19 
channels/ s ec . The input signals of the DVK were then 
digitalized and converted to tonperatures through a poly- 
nomial written according to standard tables specified by IBS 
Monograph 125 . The program for the calculator Is written 
in HP basic language and can be stored in tape cartridges. 

The accuracy of temperature measured with DAS is + 0 . 2°C . 

The other equipment used in the experiments 

includes a water bath. (COLOHA-Germcny, KLZ42.60-2-16QD) which 

o 

maintains the water temperature to an accuracy of + 0.05 C. 
Figure 5 . 2 h shows the photograph of the toot coll for two 
dimensional ca.se, water bath and DAS used for the experiments, 


5.3 EXP SRI LENT AD PROCEDURE : 

At the start of ca.ch exp orina.it al run the HP basic 
program was transferred from the tape cartridge no the 
memory of DAS. The program was written in such a way that 
the temperature lias recorded for every 15 see in case of two 
dimensional experiments and for every 1 min in caso of one 
dimensional expo rinont s . 

The horizontal and vertical locations of all the 
thermocouples wore then measured with the help of a, travelling 
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microscope ( OSAW— Indio.) to on accuracy of 0.01 mu. The 
"top surface was then fired at a known height inside the 
tost coll, The circulation of cooling rater was started 
10 min before pouring the molten war for solidification 
so that water temperature reached a constant value . The 
capillary tubes supporting the- thermocouples were heated 
to initial temperature (ih) of the molten war by blowing 
hot air on then using a hot air blower (IIOUIIIJEZ-T'rance, SII2A) . 
This was to ensure that the we:: did not solidify on the 
capillarios and the thermocouple beads, i; mediately after 
pouring it inside the test cell. 

The molten woo: heated to predetermined temperature 
T^ was then poured into the test coll through the opening 
provided at the top surface. Immediately after pouring the 
wax , the top surface* was covered with insulating mat e-rials 
and the program in DAS This run. The wax was allowed to 
solidify for 10 min in ease of smaller test cell and 2 hrs 
in ease of bigger tost cc3_l. At the end of each run, the 
unsolidified molten wax was drained out through the outlet. 

The solidified wax was then taken out for subsequent runs, 
by cutting it into pieces with a sharp knife. 


*■0 

rpi 


.ic- experiments were performed for different super- 
heats of the molten wax and for different aspect ratios of 
the test cell. 
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5.4 RESULTS; 

5,4.1 Two. D inens io nal Exp er 1 nen t s ; The 'variations of tempera— 

tare with tine measured by different themocoup les, for the 
nolten wax poured at superheats (0) of 0.52, 0.61 and 0.33 are 
shown in Figures 5*3, 5*4 and 5*5 respectively. The plots are 
drawn for only selected thermocouples out of a total of eighteen 
for clarity. The thick lines show the results fron analytical 
study described in Chapter 4. For the computations of the 
Eqs.(4*2.15) to (4.2.23) the physical properties of wax are 
chosen at a tenperature of 62 °C. Since the variation of thermal 
diffusivity a, with the temperature is very little, the constant 
property assumption is valid for comparison of the analytical 
results with the experimental results. The parameters y/Y Q and 
x/L show the horizontal and vertical positions of the thermo- 
couples fron the center and bottom of the test cell respectively. 
Fron the figures it is evident that the solidification takes 
place at a constant temperature and the slopes of the curves 
change at the solidification point. The curves lying above the 
solidification point represent the temperature distribution ii: 
the melt whereas those lying below the solidification point 
represent the tenperature distribution in the solid. Figures 
5.3 - 5.5, show that the experimental points fall reasonably 
well within the theoretical range, except during the initial 
stages of solidification. Ranachandron, Gupta and Jaluria 
(1981d) have discussed that the deviation at the initial 
stages nay be due to the turbulence created by pouring the 
nolten wax into the set up, which takes certain tine to 
subside. 
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From tho t onpcr aturo vs time plots (Figures 5.3 to 
5.5), the time at which the interface touches a particular 
tli orno couple can ho easily found out. Since the positions 
of the thermocouples are pro clot emined, the thickness of 
tho solid iomcd at different tines aro interpolated. 

Figures 5. 6-5. 8 show the interpolated experimental results 
of tho variation of solid thickness with time at different 
x/L values, and for the superheats of 0.32, 0.61 and 0.83 
respectively. The thick lines shorn the analytical results 
for comparison computed bj r Eqn. (4. 2.21d) of Chapter 4. The 
rate of solidification decreases with x/l due to the presence 
of natural convection in the molt. Tho figures show that the 
experimental points match well with the analytical solutions. 

A comparison of Figures 5.6, 5.7 and 5.8 shows that tho rate 
of solidification decreases with increasing superheat which 
is reasonable as explained in Chapters (3 and 4). The 
experimental results are compared only to a limited theoretical 
range duo to enormous computation time required (a CPU time > 

45 nin) to solve tho problem numerically for complete solidi- 
fication. 

Figures 5.9 raid 5.10 show the theoretical as well 
as the experimental temperature vs time relationship for a 
superheat of 0.66, the height L of the enclosure being 5,6cn 
(1/2 T q = 2.5, 3.0 ) respectively. Similar to the previous 
cases, the agreement between tho theory and tho experiments 
is fairly good. Figure 5.11 shows the variation of the solid 
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one can 

as 

nunc the s 

olidific; 

at ion 

to 

b u 

one dimensional if 

the asp 

eCt 

ratio is 

high. 






5-4.2 On e Dim ensional Case; The result a of the experiments 
conduct od in the- bigger tost coll for the aspect ratio of 4-7 
are compared with the theoretical ro cults conputed by one 
dimensional analysis discussed by Eanachandran, Jaluria and 
Gupta (1981a). The theoretical end orp eminent al value s of 
t enpc-ratureS measured by different thermocouples at various 
times f for the superheats of 0.18, 0.35 and 0.6 arc diown in 
Figures 5.12-5.14 respectively. Figures 5.15-5.17 show the 
thicknoss of the- solid formed at various tinu for the super- 
heats of 0,16, 0.35 and 0.6. In all these cases the agreement 
between the theory and experi: tents is quite good during the 
initial stages of solidification. As solidification progresses 
the gap between theory and experiments increases since tuo 
two dimensional effects come in. The results are not compared 
mth the two dimensional analysis as the theory developed in 
Chapter 4 is valid only for the laminar range- (la 110 6 ), 
whoros.s in case of the bigger test Cell uhe Rayleigh number 
is greater than 10 8 which lies in the turbulent zone. The 
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turbulent rn a lysis ■will include the* use of tine averaging 
terns and terns related to random fluctuations in the tlirue 
dimensions. Hence- both the analysis and the solution 
methodology trill be quite complex end involved, compared 


to the laminar analysis. 

5.5 conclusion: 

Experiments more conducted to find the- rate of 
interface movement and the t emperaturo distribution in the 

melt and i he solid during soj — LCtix.*. cation using paraffin 

*. 

wax as solidifying material. I' he* temperature of the* solidi- 
fying wax was measured imth Iron— Constantan thermocouples 
fixed inside* the test coll. The temperatures wore recorded 


using a high speed Data Acquisition System. 


The* experiments w- 


e conducted for different super- 


heats and aspect 
thickness of the 
with the theory, 
experiments was 


ratios. The* measured temperatures and tile 
solid formed at different times were compared 
The agreement between the theory and 


found to be reasonably good. 



CHAPTER 6 


CONCLUSIONS AUD RECOmiEHDATIOHS 

6.1 SUI-MAilY A m COIICIUSIOHS: 

An investigation Has been made to study the thermal 
and fluid flow effects during solidification in an enclosed 
region. The study has been divided into three sections: 

(i) Comptitations of solidification inside one 
dimensional vertical slot and vertical 
cylinder. 

(ii) Computations of two dimensional solidification 
in a rectangular enclosure with natural 
convection in the melt. 

(iii) Experimental verification of the theoretical 
results obtained vie, steps (i) and (ii). 

In the one dimensional problem, equations were 
written describing the heat transfer in the melt, the 
solid and the mold. Consideration has also been given to 
the natural convection in the melt. The momentum and 
energy equations were uncoupled and solved separately. At 
the outer surface of the mold, boundary conditions of first, 
third and fourth kinds were considered. The governing 
equations were non— ainensionalized and solved using implicit 
finite-difference technique, for various non-dimensional 
parameters such as the Prandtl number, Stefan number, Biot 


number etc 





The temperature distribution in the nelt , the 
solid and the mold, and the velocity distribution in the 
nelt were plotted. It has been found that value of the 
Grashof number is quite significant in determining the 
nagnitudo of natural convection in the nelt. In case of low 
Prandtl number fluids, appreciable fluid notion persisted 
even after the removal of ell the superheat. The magnitude 
of the fluid notion, produced by natural convection was 
found to bo quite appreciable and comparable to the terminal 
rising velocities of typical non-nctallic inclusions of 
liquid netals. The effect of various non-dimensional 
parameters on solidification has also been studied by plotting 
the rate of in terface nov orient for then. The rate of solidi- 
fication of nolten steel comput ed through present work, was 
compared with the work reported in uhe literature and the 
agreement was found to be good. 

In case of two dimensional analysis, the thermal and 
fluid motion that arise due to temperature gradients during 
solidification in a re ct angular enclosure were considered. 

The equations wore first solved without considering the 
go lidif icpition in order to have a cacck on the numerical 
scheme and later, solidification was included. Tne equations 
of continuity and momentum in the nelt were combined and 
written in terms of vorticity and stream function. These 
equations along with the energy equations in the melt, the 
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solid and the mold were no n-d in ens i o nali zed and solved 
for the following boundary conditions s 

(i) The top and bottom surfaces of the enclosure 
were kept adiabatic whereas the outer surfaces 
of the mold were kept at a constant temperature. 

(ii) The top and bottom surfaces of the enclosure 
were adiabatic, while the mold ourfcer surfaces 
exchanged energy with the surroundings b;y 
convection and radiation. 

(iii) The top surface and the mold outer surfaces 

exchanged energy by radiation and convection, 
the bottom surface being kept adiabatic. 

The vorticity and energy equations were solved 
using Alternating Direction Implicit technique. The elliptic 
stream function equation was solved using Sxiccessive Over- 
Relaxation method. The steps involved in the computation 
have been discussed in Section 4.3.6 of Chapter 4. The 
dependence of the interface movement on various non-dinensional 
parameters such as Rayleigh number, Brandt 1 number, Stefan, 
number, aspect ra/fcio, superheat, Biot number and Radiation 
constant was considered. Streamlines and isotherms were 
plotted for different Rayleigh numbers. 

It has been found that the shape of the melt region 
obtained, ruth natural convection in the molt, is markedly 
different from the pure conduction solution. When the top 
and bottom surfaces of the enclosure were kept adiabatic. 



the thickness of the nelt region varied along the length 
of the enclosure, the thickness being the smallest at 
the bottom and greatest at the top. 

fhe inclusion of natural convection in the analysis 
increased the heat transferred at the interface. This 
has been shown by plotting the ratio of the heat transferred 
at the interface with natural convection in the melt to 
the heat transferred at the intcrfa.ee without natural 
convection, at various tines. The results were computed 
for different Rayleigh numbers and compared# The conputotioi 
made to study tho effect of Prandtl number on tho interface 
motion indicated that P 3 ? >_ 10 does not have any effect 
on tho rate of movement of tho interface and an asjnaptotic 
solution was reached. A comparison of the thickness of 
the interface at different times and different aspect 
ratios indicated that the solution approaches the one 
dimensional case when the aspect ratio is high (L/2Y 2L 5). 

Experiments were performed in rectangular test 
cells using paraffin wax as the solidifying material. The 
paraffin wax was chosen duo to its low molting and easy 
handling nature. The temperature of tho solidifying wax 
was measured with Ixon-C-onst amt an thermocouples kept a.t 
predetermined levels inside the test cells and recorded 
with the help of a Hewlett Packard advanced Data Acquisition 
system. Prom the temperature vs tine plot, the position 



of tho interface was found at cliff or oat tines. 


The 


exp eriu out al resv.lt s wore conpared with the theoretical 
solutions. 

In the co.se of comparison between. the neasured and 
conputcd tenperatures, the agreement was found to be 
reasonably good ere opt during the initial stages of solidi- 
fication due to turbulence created, in pouring the woe: into 
tho test colls. The difference in thickness of the wax 
solidified, measured at different heights of the solidifying 
wax was found to be negligible in tho case of high aspect 
ratios. Similar results were obtained for tho experiments 
conducted for one dimensional case. Thus tho results 
obtained through tile experiments confirmed tho conclusion 
arrived through the two linens i oral analysis that the 
solution for the rate of solidification approaches one 
dii iensional case when the aspect ratio is high. 


5 


6 . 2 11EC QHKEJ© AT I OHS : 

In the present work, the computational program fou 
two dimensional solidification, bee:. no unstable for ila >. 10 
and Pr > 1. Modifications can be made incorporating bettor 
numerical techniques such as upwind differencing scheme, 
so that the computer program works for higher Rayleigh 
number. In the case of fla > 10°, which signifies tho natural 
convection in turbulent region ^analysis na^ bo made as 
discussed by Lin and Churchill (1S7C). This is quite essential 
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■because in case of ingot solidification of netals, the 

Hayle-igli amber calculated for a no Id of dimensions found 

12 

in typical ingot casting is of the order of 10 uhich is 
in highly developed turbulen *u n ci'u Mr? 1 convection, region. 

fjio computational m. .lysis undo for rectangular 
on closure can bo easily extruded to ti to di: icrsional solidi- 
fication pr obi o:i for cylindrical. enclosure. This can bo 
no.de uitl.out changing the :• air structure of the program 
written for rectangular enclosure o.s the only change in the 
governing equations, Eqs (4.2.5) --(4.2.10), is the additional 
convection terns due to ti a introduction of r- coord in ate 
as slioiei below; 


Energy equation in nolt? 

92 1 9T 1 91 'l ,1 ' r ^1 W 

~ri + u *“ i + v Tr = a i r — ' Cr “n? } + 


3X 

( 6 . 2 . 1 ) 


Vorticity: 


s “ 3 u 3 tu - a) 2 u \ n jg- i 

4- v - V ~’- J = v ( V 0) - o)~SP - 1 - 

r r r 


Tt + w ~ 


v •- 


3^ 

( 6 . 2 . 2 ) 

where V 2 is the Laplace operator in cylindrical coordinates, 


St r can functions 


a 4 


0) 


i ( _iX + i .-11 + 

r ^ 3r 2 r sr ^ 


(6.2.3) 



141 


Velocity: 

- 2 1 - 1 
U ~ r 3r’ ~ r "ITS 

Energy equation in solid: 


( 6 ; 2i4) 


s /-l a / 

Tt = a s ( r “? r 


■§) (6.2.5) 


Energy equation in no Id! 


a ri^ r 


Jj) + 

T J 


,2 1 
? n\ 


( 6 . 2 . 6 ) 


Tlio analysis for the two d incus ional solidification 
nay loo further extended to alloy cystous. This would find a 
direct application in nost of the netallurgical industries, 
where the nelt solidifies over a range of tonperature. In 
case of alloys, separate non certain and energy equations will 
have to bo written for the mushy zone apart from the existing 
equations for the nelt, the so3.id and the nold. This part 
can a.lso be introduced in the existing computer program 
without any na.jor change in the program structure. 

Finally, the analysis can bo extended to throe 
din elision all solidification. The introduction of the third 
dimension coupiica/b co the solution methodology and special 
methods are to be looked into to solve the three dimensional 
p a r ab o li c e qu at ions. 

In the case of experiments, apart from the t enp eratur e 
neasurononts, the flow developed in the nelt can be studied 



with ol ic help of tracer particles end optical nothods as 


discussed by Szokels' end Jaosal (1978). This is possible 
only if the solidification occurs below room t cop creature, 
so that the material does not solidify on the Perspex walls 
of the test cell, immediately after pouring the licit and 
thereby obstructing tin view for any optical noasurenent . 
For this purpose any low ielting n~paraffin between C-^H^q 
raid C^gH^o no -7 ho need. The flow developed in the no It 
can be observed using Sclilioren photographic technique end 
the interface movement era be phot ographed . For experiment 
to study the solidification of binary alloy systems, a 
solution of ammonium chloride in water con bo used as the 


solidifying material with methanol 


circulated in the 


jacket 


of the test coll as coolant. Experiments ray also be 
conducted, for other boundary conditions described in 
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APPEitDIX I 


SIKPLLET .CATI 017 OP I10ME1ITUM EQUATION 


The momentum equation Tor the u—componeiii is 


vrxx'cen as 


^ .jlLb 

" »y 2 


Og ~ 


( 1 . 1 ) 


The pressure term appearing in the PJIS resy be 
broken into tiro toms, one aue to hydrostatic pressure p 

Ci 

i ’ St 

and the other -duo to motion of fluid. The variation of 
along x-dire cti on is neglected, be cans e of the one dimensional 
assumption. Since the gravitational force is acting on the 
x-dire ct ion, the hydrostatic pressure term p_ can be written 


p sai s 


( 1 . 2 ) 


Sub st itut ing (1.2) in ( 1 . 1 ) 


D 3 u * u 

p Tt = " 2 

^ y 


+ g ( p + - p ) (1.3) 


The densities appearing in Ban. (1.3) can be related 
to the temperature through the Bouscinesq. approximation, 
discussed by Jaluria (1980). 

According to Boussinosq. approximation, density 
variations are considered only insofar as they contribute 
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to buoj'ancy , out otherwise neglected, 
average temperature of the fluid and 
ponding density, then, 


If T is the 
P is the corres' 

Sal 


sat 


-PM' 


T) 


(1.4) 


where p = t sniper atu.ro coefficient of cubic expansion. 


i'hus t Lo Eqn . ( 1 . 1 ) be- cores. 


3 u 

Vc 


v 


— ^ + gp (I 
3 7 


(1.5) 



APPE1JDIX II 


RADIATIVE B QUID ARY G 0111 HI Oil 

Prom the Eqn. (3.2.151) , the radiative boundary 
cond. it i o n is vri tt on a s , 

n** 


3 T 


m 


3 y * 


- lio l [ W* + Ip 4 - O 


(II.1) 


This con be simrlif iod to 




n Oi(^ 2 + 1 o 2 ) K + (“-2) 


where 


T, 


# 


*¥■ 

T + T 
c m 


In finite-difference- form , the above equation nay be 
written as. 


y. 

n , ;i , n-f 1 m , ,i — i ,n+l __m* j L c + l*i(T* . +T*) 

T^T - j,nfl XuC l U -t,3,n + x c ;u t, 3 ,n c ; 


m 


(II. 3) 


Since T* . is kiiomi from the previous tunic step, Eqn* (11*3) 

t f 3 * 33. 

can "be simplified to i 


(1 + A Q > “ ■ i 'n,j-l,n+l 


0 


(II. 4) 


where £ c - ac 1 + ?') ("■&,;) ,n + 

Ean.(11.4) 4 s linear and con be solved along with the other 
oo nations formed by implicit, finite-difference method. 
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AP PUD IX III 

PP-OPB.li'X ES OF STEEL 

Pron Elliott and Gleisor (I960), the following 
properties are t alee ns 

Specific heat 0 = 6.5 r 10^ J/K q °K 

Jr 

Density p = 7.1 r. 10^ Z<s/n? 

—*7 P 

iCinemaxic viscosity v = 9.4 r 10“ n / s 

Prandtl nunber Pr = 0.15 

latent heat ! 0 = 226 x 10^ Sj Kg 

Iherual conductivity k-^ = 0.464 J/s m °K 

6 2 

Thermal diffusivity a 1 = 6.27 n 10" n/s 

Coefficient of cubical 

eurpansion (3 = 1.5 n 10 ^ °K ^ 
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APPEKDIX IY 


IIOYIiTG B GUI DAIlY ECU AT I Oil 


Tho energy balance at the melt-solid interface is 


written as, 


k. 


?T 0 SI. 

' an - k l = P L q V 11 


(IV. 1) 


whore n = outward normal to the interface (i,o. into the 
melt ) 

v = velocity of tho interface in the normal direction 

Since tho interface is an isotherm tho following condition 
is written, 


T 


1 


'sat 




s 


(IY. 2) 


The interface position can ho defined by the following 
equation 

f’(r,y,z;,t) = 0 (IV. 5) 

In tho case of two dimensional solidif ication only 
the effect of x and y directions is considered. liorcovor, 
the temperature gradients at tho interface in tho r-diroct ion 
can Too expressed in terns of gradients in y-direction by the 
differentiation of Eqs (IV. 2) and (IV. 5) as 


9f 


_1 
a x 


a-pt 


3f*/ay 


9 T. 


ay 


(IV, 4) 
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and. 


91 


?T 


_s _ 3f l / s 

is ~ a£!/1y "a*F 


(IV. 5) 


Equation (IY.1) can le expressed in the following 
form as discussed "by Patel (1968). 


k D ,S o vf -k x ^v f .=- (W.6) 

Substitution of Los. (IV. 4) and (IV. 5) into Sqn.(IV,6) 
results in tlio fo Hearing crprcssioii 




a l 


9 1. 


__H _ k 
ay i 


a y 


-] 


T af Vat 

- P \ -nyVy 


(IV. 7) 

If f»(x, 3 ',t) = y-Y T (:c,t) = 0, Eon. (IV. 7) reduces to 


<*□ TF 


* _a )( i +( 8 r.) 2 ) 

1 a y A ^ ar. ' ' 


= PL 


a i 1 

q at . 


(IV. 0) 




APPIdID] 


V 


c oo.rj) D\as 2 TiiAiisroidi&a ion 


Tlie original variables r, y and t cm be replaced 


by the 


u .i. r alb 


s formed variables 


* * 


y and r 


ip plying tlie 


non-dinen si on al parar.eteis described in Chapter 4 (Eqn. 
(4-2.14)). 

Thus , 


For nelt : 


3 

Tt 




^or solids 


3 y 


a 

Tt 


a _i 

22 
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1 
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1_ 

T* 


* 

r l 


n 




i 1 


Zi 

Y* 


3 y< 

at 7 


a 22 

3:-:. 


3y 


1 


3 7 
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* 

'1 


3 Y* 


on 3 y* - 1 

4 — * + # 

U* -o X . ?y* 


(Y. la) 

(Y.lb) 

(V.lc) 


(Y. 2a) 
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l 
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r+ ( 


v s 


* 1 \J 


■,) 


3 y 


U i -Y ,; 
o 


3 y, 


3 Y< 


3 y 


(V. 2b) 
(V . 2c) 


Applying Eqc (Y-la) - (V-2c) to Eqs (4. 2.5 )-(4. 2.10) of 
Chapter 4, the following transformed equations are obtained 
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The Eqs( 7 . 2 )-(V.G) ere quite coupler and difficult 
to solve even with a numerical scheme as the convergence 
of the iterative scher:e is quite liheljr to he affected* 

In the actual iornulation it has teen assunec. that the 
interface remains stationary for the period in which heat is 
transferred across it Jhuring this period, the equations 
describing the i .omentum and energy of the melt and the energy 
equations of the solid and mold are solved assuring pseudo- 
st eady state as discussed at great length by Sparrow , 

Patankar and Ranadliyani (3.977). Because of tnis assumption^ 
the tern containing 3 6 / Ot* in the -CHS or Pqs (7.3) * (7.4) 
and (Y.7) can be neglected. Further, it has also been 
assumed that the thickness of the melt region 6 varies 
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slowly with, x* so that the tes:.& in vo Ivin 
2 2 

a 6/ a::* appearing in '■'lie 1UIS of Eqs.(V 
neglected . Sparrow, Pat an far and haiiaehy 
discussed that in the case of two dir.enci 
this assumption is quite valid as a first 
nost cases of phase change problems with 
in the melt . 

3 e caus e o f t h e e e as eurrpt ions Eos . ( 
simplified to Eqs. (4.2.15) to (4.2.20) £ 


■Z S 6/ 3 n* and 
.3) to (Y.7) are 
uni (1977) have 
o;ial so lidif t a sb ion , 
apurosi nation in 
natural, conv e ct ion 

Y.3) - (V.8) are 
:iven in Chapter 4- 
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SBCQffl) TJPVIITD D ZPTEPJ].:TC III G - I-lhiHOD 
The vorticity equation is written as. 
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Tlie second upwind differencing scheme is -used to 

■write the convective terms appearing in HIS of Eqn.(VX.l) 
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in finite-difference form, 'finis, for example, U* is 
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In. the first half of the tine step, Eqn.(YIIl) is written aSj 
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This con he simplified to, 
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where El = ~ (“~^) (^f ) 2 

E 2 = (.JLL. w yt) 

- Ay 

E 3 = 2 E 1 

E4 = (lLkl i l)(^ )( |_) 
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Similarly for the second half of the tine stop, Eqn.{VII«l) 
can he cimpliiied to; 
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The subscripts 1 , 2,3 appearing in Eqns, (VII* 2 ) ~ (YII* 4 ) 




are the variables 
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APP ilTD IZ V IIi 
PPOPBRTIB S 01' P AllAPFIIT WAX 

Prom Bennett (1559) and kelson (1958), fio following 
properties are taken % 
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